Bayesian data analysis of programming languages

Part I: FSE data

1 Data preparation and descriptive statistics

We use the original data from the FSE publication (Ray et al. 2014).

First, manually do \(\log\) transformation,

fse.data$commits_log <- log(fse.data$commits)
fse.data$insertions_log <- log(fse.data$insertions)
fse.data$max_commit_age_log <- log(fse.data$max_commit_age)
fse.data$devs_log <- log(fse.data$devs)

Then, remove columns not needed and add a column that represents the ‘project_id’ as a numeric.

fse.data = subset(fse.data, select = -c(domain))
fse.data$project_id <- as.integer(fse.data$project)

Our data frame now looks like this,

glimpse(fse.data)

## Rows: 1,127
## Columns: 13
## Groups: project [729]
## $ project            <fct> _s, 4clojure, 4clojure, accelerate, ack, ActionBar…
## $ language           <fct> Php, Clojure, Javascript, Haskell, Perl, Java, Jav…
## $ commits            <int> 174, 629, 75, 985, 97, 213, 796, 39, 64, 1404, 42,…
## $ insertions         <int> 3232, 9180, 39240, 99003, 425, 14657, 225043, 1628…
## $ max_commit_age     <int> 696, 774, 509, 1644, 2031, 192, 874, 1358, 1297, 1…
## $ n_bugs             <int> 0, 128, 29, 165, 12, 68, 157, 11, 9, 397, 19, 71, …
## $ devs               <int> 40, 21, 10, 17, 6, 5, 42, 13, 9, 170, 6, 35, 24, 1…
## $ language_id        <dbl> 13, 4, 10, 8, 12, 9, 9, 5, 10, 15, 13, 8, 11, 11, …
## $ commits_log        <dbl> 5.159055, 6.444131, 4.317488, 6.892642, 4.574711, …
## $ insertions_log     <dbl> 8.080856, 9.124782, 10.577452, 11.502905, 6.052089…
## $ max_commit_age_log <dbl> 6.545350, 6.651572, 6.232448, 7.404888, 7.616284, …
## $ devs_log           <dbl> 3.6888795, 3.0445224, 2.3025851, 2.8332133, 1.7917…
## $ project_id         <int> 1, 2, 2, 3, 4, 5, 6, 7, 7, 7, 8, 9, 10, 11, 12, 13…

where ‘n_bugs’ (\(\mathbb{N}^+\)) is our outcome variable and ‘language_id’ (\(\mathbb{N}^+\)), ‘project_id’ (\(\mathbb{N}^+\)), ‘commits_log’ (\(\mathbb{R}^+\)), ‘max_commit_age_log’ (\(\mathbb{R}^+\)), ‘devs_log’ (\(\mathbb{N}^+\)), and ‘insertions_log’ (\(\mathbb{R}^+\)), are our potential predictors.

We should check to see that we have no NAs nor zero-inflation in the outcome variable:

table(is.na(fse.data))

## 
## FALSE 
## 14651

table(fse.data$n_bugs == 0)

## 
## FALSE  TRUE 
##  1125     2

Since n_bugs \(\in \mathbb{N}^+\) we would expect to use a Poisson(\(\lambda\)) likelihood (and of course we tried that first). If we, however, look at some more descriptive statistics of ‘n_bugs’ we see that there is a large difference between the mean and the variance,

summary(fse.data$n_bugs)

##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
##      0.0     20.0     65.0    500.7    228.0 121302.0

var(fse.data$n_bugs)

## [1] 15031006

clearly indicating that we need to model the variance separately (each Poisson count observation should have its own rate). Hence, we’ll assume that the underlying data-generation process approximately follows a negative binomial (Gamma-Poisson) distribution \(\mathrm{NB}(\lambda,\phi)\).¹

2 Initial model design

It is prudent to start model design by building a small model and then move to more complex models. But what do we know so far? Well, the outcome (dependent variable), \(y\), is a count starting from \(0\) (albeit only with \(2\) zeros). We also know that the mean is considerably smaller than the variance indicating that the underlying data generative process is a Negative-Binomial.

As predictors (independent variables) we have project numbers (unique ID for each project), the language used (language_id) and then transformed (\(\log\)) predictors: maximum commit age (max_commit_age), number of developers (devs), number of insertions (insertions), and number of commits (commits).

2.1 Take 0

Very often, when designing statistical models, we want to have a simple model, \(\mathcal{M}_0\), to compare against. A simple model in this case could be a model where we simply estimate the mean for each language \(\alpha_{\mathrm{LANG}[i]}\):

\[\begin{equation} \begin{aligned} \mathrm{n\_bugs}_i & \sim \textrm{Negative-Binomial}(\lambda_i,\phi) \\ \log(\lambda_i) & = \alpha_{\mathrm{LANG}[i]} \\ \alpha_j & \sim \mathrm{Normal}(\bar{\alpha},\sigma_l) & \textrm{for LANG}=1,..,17\\ \bar{\alpha} & \sim \mathrm{Normal}(0,5)\\ \sigma_l & \sim \mathrm{Weibull}(2,1)\\ \phi & \sim \mathrm{Gamma}(0.01,0.01) \end{aligned} \end{equation}\]

There’s nothing strange with the above. We simply assume that n_bugs is distributed according to a Negative-Binomial process with the rate \(\lambda\) (just like in Poisson) and a parameter \(\phi\) that controls the variance (which must be positive). For the linear part in our model (Line 2) we use a \(\log\) link to translate from the probability space \(p \in (0,1)\). Let’s look what the priors imply.

2.1.1 Prior predictive checks

p <- get_prior(n_bugs ~ 1 + (1 | language_id),
           family = negbinomial(),
               data=fse.data)

p$prior[1] <- "normal(0,5)"
p$prior[2] <- "exponential(1)"

m0p <- brm(n_bugs ~ 1 + (1 | language_id),
           family = negbinomial(),
               data=fse.data,
               prior = p,
               sample_prior = "only"
               )

y <- fse.data$n_bugs # our empirical outcomes
y_rep <- posterior_predict(m0p) # samples from the posterior 

ppc_dens_overlay(y, y_rep[1:100,]) + 
  scale_x_log10(
   breaks = scales::trans_breaks("log10", function(x) 10^x),
   labels = scales::trans_format("log10", scales::math_format(10^.x))
 ) +
  annotation_logticks(side = "b", outside=TRUE) +
  coord_cartesian(clip = "off") + 
  xlab("Number of bugs") + 
  ylab("Density") +
  geom_vline(xintercept = max(fse.data$n_bugs), linetype = "longdash")

The prior predictive check shows that the priors show great variability (the \(x\)-axis is log10() transformed), and we also have a number of values that are \(+\infty\). The dark line is our empirical data \(y\), while draws from our posterior (which has only been sampled using the priors) are the light blue lines (\(y_{\mathrm{rep}}\)). The dashed line indicates the maximum number of bugs in our empirical data \(y\). But what maximum values (that is not \(+\infty\)) do the priors provide?

max(y_rep, na.rm = T)

## [1] 1161985498

That’s a very large value.

2.1.2 Posterior predictive checks

Running the above model with default priors and settings (in later examples we will set our own priors).

m0 <- brm(n_bugs ~ 1 + (1 | language_id),
          family = negbinomial(),
          data=fse.data,
          prior = p,
          iter = 5e3,
          refresh = 0
          )

We then conduct posterior predictive checks.

pp_check(m0) + 
  scale_x_log10(
   breaks = scales::trans_breaks("log10", function(x) 10^x),
   labels = scales::trans_format("log10", scales::math_format(10^.x))
 ) +
  annotation_logticks(side = "b", outside=TRUE) +
  coord_cartesian(clip = "off") + 
  xlab("Number of bugs") + 
  ylab("Density")

As we can see (note that the \(x\)-axis is logarithmic base \(10\), with log ticks displayed along the axis), the model struggles and we don’t have a good fit. We’ll add the other variables at our disposal and see if we can improve things. But first, let’s have a look at the diagnostics.

2.1.3 Diagnostics

A divergent transition indicates that the Hamiltonian trajectory has departed from the ‘true’ trajectory and we don’t want this to happen.²

np <- nuts_params(m0)
sum(subset(np, Parameter == "divergent__")$Value)

## [1] 0

The \(\widehat{R}\) is the ratio of within- and between-chain variance. As a general rule of thumb, a stationary posterior probability distribution is indicated by \(\widehat{R} < 1.01\), i.e., it should approach \(1\) when \(N \rightarrow \infty\). We see that the diagnostics here is just below \(1.01\).

# rm NAs since they have not been estimated
max(rhat(m0), na.rm = TRUE)

## [1] 1.001869

Next, we should check the effective sample size (ESS) for each parameter we estimated. Generally speaking, we always want \(\mathrm{ESS} \gtrapprox 0.1\), i.e., \(10\)% of the total sample size at a minimum, for each parameter we estimate (or at least it should be in the hundreds in absolute values).

# rm NAs since they have not been estimated
min(neff_ratio(m0), na.rm = TRUE)

## [1] 0.1126871

Finally, we visually inspect trace plots to see that they look like hairy caterpillars (i.e., that they mix well),

mcmc_trace(m0) + 
  theme(text = element_text(size=6))

2.2 Take I

In the original study the authors used four population-level effects: devs_log, max_commit_age_log, commits_log, and insertions_log. If we would follow a principled Bayesian approach (and for that matter a frequentist approach) many would claim that some type of variable selection should take place (since we are after good predictions). Later down below, Sect. 5, we will show how this can be done, but we do not include it here, since our intention is to follow the original study.

Our next model, \(\mathcal{M}_1\), is a model that estimates a grand mean, \(\alpha\), four population-level effects (\(\beta_1,\ldots,\beta_4\)), and then employ varying intercepts according to language.

In math speak this is what it looks like without the priors,

\[\begin{equation} \begin{aligned} \mathrm{n\_bugs}_i & \sim \textrm{Negative-Binomial}(\lambda_i,\phi)\\ \log(\lambda_i) & = \beta_1 \cdot \mathrm{devs\_log}_i + \beta_2 \cdot \mathrm{max\_commit\_age\_log}_i\\ & + \beta_3 \cdot \mathrm{commits\_log}_i + \beta_4 \cdot \mathrm{insertions\_log}_i\\ & + \alpha_{\mathrm{LANG}[i]} \\ \end{aligned} \end{equation}\]

2.2.1 Prior predictive checks

First, we’ll make sure to check what priors such a model expects us to set,

(p <- get_prior(n_bugs ~ 1 + devs_log + max_commit_age_log + commits_log + 
                  insertions_log + (1 | language_id), 
                data = fse.data, 
                family = negbinomial
                )
 )

For our \(\beta\)s we have unsuitable uniform priors, \(\mathrm{Uniform}(-\infty,\infty)\). For the intercepts, \(\alpha\), a Student’s \(t\) distribution with \(\nu=3\), and the same for the standard deviation, \(\sigma\), used for the varying intercepts (however, in this case it’s automatically truncated with a lower bound of \(0\) since \(\sigma\) cannot be negative). Finally, there is a default prior for the shape, \(\phi\), i.e., \(\gamma(0.01,0.01)\).

It might be prudent to set our own priors to see how the priors look like on the outcome scale.³ We’ll set \(\mathrm{Normal}(0,5)\) for our intercepts, \(\mathrm{Weibull}(2,1)\) for our dispersion \(\sigma\), and \(\mathrm{Normal}(0,0.5)\) for each of our four \(\beta\)s. In math notation we now have this model,

\[\begin{equation} \begin{aligned} \mathrm{n\_bugs}_i & \sim \textrm{Negative-Binomial}(\lambda_i,\phi) &\\ \log(\lambda_i) & = \beta_1 \cdot \mathrm{devs\_log}_i + \beta_2 \cdot \mathrm{max\_commit\_age\_log}_i &\\ & + \beta_3 \cdot \mathrm{commits\_log}_i + \beta_4 \cdot \mathrm{insertions\_log}_i &\\ & + \alpha_{\mathrm{LANG}[i]} &\\ \alpha_j & \sim \textrm{Normal}(\bar{\alpha},\sigma_l) & \textrm{for LANG}=1,..,17\\ \bar{\alpha} & \sim \mathrm{Normal}(0,5) \\ \sigma_l & \sim \mathrm{Weibull}{(2,1)} \\ \beta_1,\ldots,\beta_4 & \sim \textrm{Normal}(0,0.5)\\ \phi & \sim \textrm{Gamma}(0.01,0.01) & \end{aligned} \end{equation}\]

The \(\mathrm{Gamma}\) prior on \(\phi\) allows only positive real numbers, as does the \(\mathrm{Weibull}(2,1)\) for our varying intercept. Additionally, having \(\mathrm{Normal}(0,0.05)\) on our \(\beta\)s might seem tight, but the additive terms imply a variance of \(\exp^{(0.05*4)^2}=1.04\). Let’s see what the prior predictive checks will show.

Sample the model using no data.

p$prior[1] <- "normal(0,0.5)"
p$prior[6] <- "normal(0,5)"
p$prior[7] <- "weibull(2,1)"

m1p <- brm(n_bugs ~ 1 + devs_log + max_commit_age_log + commits_log +
                 insertions_log + (1 | language_id),
          family = negbinomial(), 
          prior = p, 
          data=fse.data, 
          chains = 1, 
          iter = 1000, 
          sample_prior = "only", 
          refresh = 0
          )

We used only one chain and we sampled for only \(1000\) iterations.⁴ Half of the iterations are thrown away (warm-up), which leads us to having \(500\) samples in our prior probability distribution.⁵

Let’s now plot samples, \(y_{\mathrm{rep}}\), from our prior probability distribution, and compare them to our empirical data, \(y\).

y_rep <- posterior_predict(m1p)

ppc_dens_overlay(y, y_rep[1:100,]) + 
  scale_x_log10(
   breaks = scales::trans_breaks("log10", function(x) 10^x),
   labels = scales::trans_format("log10", scales::math_format(10^.x))
 ) +
  annotation_logticks(side = "b", outside=TRUE) +
  coord_cartesian(clip = "off") + 
  xlab("Number of bugs") + 
  ylab("Density") +
  geom_vline(xintercept = max(fse.data$n_bugs), linetype = "longdash")

In short, the priors are all over the place (and a number of them are \(\infty\)). The maximum number of bugs in our dataset is depicted by the dashed line, and the maximum number our priors provide (i.e., \(\neq \infty\)) is,

max(y_rep, na.rm = T)

## [1] 265676267

Let’s sample the model with our empirical data by simply removing the statement sample_prior = "only".

m1 <- brm(n_bugs ~ 1 + devs_log + max_commit_age_log + commits_log +
                 insertions_log + (1 | language_id), 
          data = fse.data, 
          family = negbinomial,
          prior = p, 
          refresh = 0
          )

2.2.2 Posterior predictive checks

pp_check(m1) + 
  scale_x_log10(
   breaks = scales::trans_breaks("log10", function(x) 10^x),
   labels = scales::trans_format("log10", scales::math_format(10^.x))
 ) +
  annotation_logticks(side = "b", outside=TRUE) +
  coord_cartesian(clip = "off") + 
  xlab("Number of bugs") + 
  ylab("Density")

Clearly the model fits the empirical data better. This is an indication that adding the four population-level parameters makes a difference.

2.2.3 Diagnostics

Let’s check a number of diagnostic: posterior predictive checks, divergences, \(\widehat{R}\), effective sample size (ESS), and trace plots.⁶

Divergent transitions:

np <- nuts_params(m1)
sum(subset(np, Parameter == "divergent__")$Value)

## [1] 0

\(\widehat{R}\) ratio:

max(rhat(m1), na.rm = TRUE)

## [1] 1.007007

Effective sample size (ESS):

min(neff_ratio(m1), na.rm = TRUE)

## [1] 0.1922557

Traceplots of all parameters:

mcmc_trace(m1) + 
  theme(text = element_text(size=6))

Once the model has passed all diagnostics we can start putting some trust in it. Let’s see what the model has estimated!

2.2.4 The fallacy of placing confidence in point estimates

summary(m1)

##  Family: negbinomial 
##   Links: mu = log; shape = identity 
## Formula: n_bugs ~ 1 + devs_log + max_commit_age_log + commits_log + insertions_log + (1 | language_id) 
##    Data: fse.data (Number of observations: 1127) 
## Samples: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup samples = 4000
## 
## Group-Level Effects: 
## ~language_id (Number of levels: 17) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.20      0.04     0.14     0.30 1.01      931     1698
## 
## Population-Level Effects: 
##                    Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## Intercept             -1.92      0.12    -2.16    -1.68 1.00     2029     2917
## devs_log               0.05      0.01     0.03     0.08 1.00     3133     3144
## max_commit_age_log     0.04      0.02     0.01     0.08 1.00     3416     2973
## commits_log            0.99      0.01     0.96     1.01 1.00     2387     2649
## insertions_log         0.02      0.01     0.00     0.04 1.00     2520     2715
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## shape     7.75      0.39     6.99     8.53 1.00     3806     2884
## 
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Above we see the estimated \(\sigma\) for the group-level effects, i.e., the \(17\) languages which we used as varying intercepts. Then we see our \(\alpha\) (Intercept), which is a population-level effect (fixed effect), just like our four \(\beta\) estimates. Finally, we have the shape estimated. Except for the shape, all other parameters must be transformed using the exponential, since we employed a \(\log\) link function in our model to translate to the probability space.

The shape parameter, \(\phi\), is an indication of the level of aggregation/clustering among our languages, i.e., when it goes towards infinity it corresponds to an absence of aggregation among our languages. It’s also an indirect measure of variance, since \(\sigma^2 = \lambda + \lambda^2 / \phi\).

In order for us to get some estimates on the different languages we need to look into our varying intercepts, i.e., the languages. However, simply looking at point estimates to check if they are significant or not, is quite simply not sane when we have a posterior probability distribution.

ranef(m1)

## $language_id
## , , Intercept
## 
##       Estimate  Est.Error         Q2.5       Q97.5
## 1   0.18082141 0.06639260  0.050364479  0.31228134
## 2   0.04897556 0.07249342 -0.093952804  0.18717831
## 3   0.23940222 0.06559268  0.114715949  0.37018122
## 4  -0.24929348 0.07443846 -0.394913675 -0.10892506
## 5   0.08345173 0.07382626 -0.059232665  0.23097126
## 6   0.02502006 0.07390779 -0.122409477  0.17439957
## 7  -0.06819667 0.08096160 -0.224857915  0.08953470
## 8  -0.19772370 0.07749466 -0.353692659 -0.04759631
## 9   0.01359081 0.06802235 -0.117230200  0.14971426
## 10  0.09924367 0.05965052 -0.014500928  0.21802844
## 11  0.18982778 0.07174194  0.051082856  0.33527950
## 12 -0.05903588 0.09089232 -0.231973213  0.12568932
## 13  0.15698436 0.07108923  0.017180808  0.29648176
## 14  0.12828689 0.06426550  0.004135146  0.25672383
## 15 -0.08359171 0.06976338 -0.219147029  0.05362552
## 16 -0.17513121 0.07198830 -0.320555675 -0.03657739
## 17 -0.34976585 0.07716196 -0.504748351 -0.20040922

As we can see above, quite many of our \(17\) languages are ‘significant.’ But, what does it mean in practice? Well, we have a posterior probability distribution we can ask questions to. As an example, assume that we use the empirical mean for the four predictors we’ve set as population-level effects, i.e., devs, max_commit_age, insertions, and commits, then we simply check what the difference is between two of the languages, let’s take the first two, which are C and C#.

nd <- data.frame(
  language_id = seq(1:2),
  insertions_log = log(mean(fse.data$insertions)),
  commits_log = log(mean(fse.data$commits)),
  max_commit_age_log = log(mean(fse.data$max_commit_age)),
  devs_log = log(mean(fse.data$devs))
)

pp <- posterior_predict(m1, newdata = nd)
summary(pp[,1]) # C

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    81.0   365.8   470.0   494.0   598.0  1387.0

summary(pp[,2]) # C#

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    64.0   311.0   404.5   423.6   514.0  1287.0

There’s clearly a difference in number of bugs between the C (upper) and the C# (lower) summaries. A difference we will find among many languages and which we’ll explore later.

2.3 Take II

In the previous section we created a model, \(\mathcal{M}_1\), that we can continue to build on. The logical next step is to keep our varying intercepts, but add the concept of varying slopes. In short, not only will each language’s intercept be modeled separately (i.e., where it crosses the \(y\)-axis), we will also model each language’s slope separately (i.e., the ratio of the change along the \(y\)-axis to the change along the \(x\)-axis). Mathematically speaking we thus have,

\[\begin{equation} \begin{aligned} \mathrm{n\_bugs}_i & \sim \textrm{Negative-Binomial}(\lambda_i,\phi)\\ \log(\lambda_i) & = \alpha_{\mathrm{LANG}[i]} + \beta_{\mathrm{LANG}[i]} \cdot \mathrm{devs\_log}_i + \gamma_{\mathrm{LANG}[i]} \cdot \mathrm{max\_commit\_age\_log}_i \\ & + \delta_{\mathrm{LANG}[i]} \cdot \mathrm{commits\_log}_i + \epsilon_{\mathrm{LANG}[i]} \cdot \mathrm{insertions\_log}_i + \zeta_{\mathrm{PROJ}[i]}\\ \begin{bmatrix} \alpha_{\mathrm{LANG}} \\ \beta_{\mathrm{LANG}} \\ \gamma_{\mathrm{LANG}} \\ \delta_{\mathrm{LANG}} \\ \epsilon_{\mathrm{LANG}} \end{bmatrix} & \sim \mathrm{MVNormal} \left( \begin{bmatrix} \alpha\\ \beta\\ \gamma\\ \delta\\ \epsilon \end{bmatrix} , \Sigma \right)\\ \Sigma & = \left( \begin{matrix} \sigma_\alpha & 0 & 0 & 0 & 0\\ 0 & \sigma_\beta & 0 & 0 & 0 \\ 0 & 0 & \sigma_\gamma & 0 & 0 \\ 0 & 0 & 0 & \sigma_\delta & 0 \\ 0 & 0 & 0 & 0 & \sigma_\epsilon \end{matrix} \right) \mathbf{R} \left( \begin{matrix} \sigma_\alpha & 0 & 0 & 0 & 0\\ 0 & \sigma_\beta & 0 & 0 & 0 \\ 0 & 0 & \sigma_\gamma & 0 & 0 \\ 0 & 0 & 0 & \sigma_\delta & 0 \\ 0 & 0 & 0 & 0 & \sigma_\epsilon \end{matrix} \right) \\ \alpha,\ldots,\epsilon & \sim \mathrm{Normal}(0,0.5)\\ \sigma_\alpha,\ldots,\sigma_\epsilon & \sim \mathrm{Weibull}(2,1)\\ \mathbf{R} & \sim \mathrm{LKJ}(2)\\ \zeta_j & \sim \mathrm{Normal(\bar{\alpha}, \sigma_p)} & \textrm{for PROJ} = 1,..,729\\ \bar{\alpha} & \sim \mathrm{Normal(0,5)}\\ \sigma_p & \sim \mathrm{Weibull(2,1)}\\ \phi & \sim \textrm{Gamma}(0.01,0.01) \end{aligned} \end{equation}\]

There are a few things that makes this model look different compared to the previous model. The first thing is that for each population level effect we set a varying slopes (Lines \(2\)–\(3\)). Then, on Lines \(4\)–\(5\) we model the slopes using a multivariate Normal (we want to model the correlation between slopes). We also set a Lewandowski-Kurowicka-Joe (LKJ) prior to model the covariance matrix of the multivariate normal distribution (Line \(9\)). Finally, as seen at the end of Line \(3\), we add another varying intercept (\(\zeta_{\mathrm{PROJ}[i]}\)) according to project ID (and on Line \(10\) we set the prior with hyperparameters).

The idea with the latter part is that we believe projects could learn from each other by employing partial pooling. Some projects have several rows of data, while some only a few. If the projects having more data informs projects with less data we might learn more, but still avoid overfitting.

Before we continue with our prior predictive checks it might be worthwhile to see what an \(\mathrm{LKJ}()\) prior implies. This prior is used for covariance matrices exclusively. The difference between \(\mathrm{LKJ}(1)\) and \(\mathrm{LKJ}(2)\) is that in the former case it implies a virtually flat prior, while in the latter case we claim that we are skeptical of extreme correlations. In our case, we have no better information than that.⁷

Figure 2.1: Correlation priors LKJ(1) (black) and LKJ(2) (blue)

2.3.1 Prior predictive checks

As usual, we’ll first check the default priors and then set our priors as stipulated in the previous section.

(p <- get_prior(n_bugs ~ 1 + devs_log + max_commit_age_log + commits_log + insertions_log +
              (1 + devs_log + max_commit_age_log + commits_log + insertions_log | language_id) +
              (1 | project_id),
          family = negbinomial,
          data = fse.data
          )
 )

p$prior[1] <- "normal(0,0.5)"
p$prior[6] <- "lkj(2)"
p$prior[8] <- "normal(0,5)"
p$prior[9] <- "weibull(2,1)"

Next we sample our model using only the priors and conduct prior predictive checks.

m2p <- brm(n_bugs ~ 1 + devs_log + max_commit_age_log + commits_log + insertions_log +
            (1 + devs_log + max_commit_age_log + commits_log + insertions_log | language_id) +
            (1 | project_id),
          family = negbinomial(), 
          prior = p, 
          data=fse.data, 
          chains = 1, 
          iter = 1000, 
          sample_prior = "only", 
          refresh = 0
          )

y_rep <- posterior_predict(m1p)

ppc_dens_overlay(y, y_rep[1:100,]) + 
  scale_x_log10(
   breaks = scales::trans_breaks("log10", function(x) 10^x),
   labels = scales::trans_format("log10", scales::math_format(10^.x))
 ) +
  annotation_logticks(side = "b", outside=TRUE) +
  coord_cartesian(clip = "off") + 
  xlab("Number of bugs") + 
  ylab("Density") +
  geom_vline(xintercept = max(fse.data$n_bugs), linetype = "longdash")

So, looks good. The priors are all over the place and additionally allow quite absurd values. Let’s now sample the model with data.

m2 <- brm(n_bugs ~ 1 + devs_log + max_commit_age_log + commits_log + insertions_log +
            (1 + devs_log + max_commit_age_log + commits_log + insertions_log | language_id) +
            (1 | project_id),
          family = negbinomial(), 
          prior = p, 
          data=fse.data, 
          refresh = 0, 
          control = list(adapt_delta=0.95)
          )

2.3.2 Posterior predictive checks

pp_check(m2) + 
  scale_x_log10(
   breaks = scales::trans_breaks("log10", function(x) 10^x),
   labels = scales::trans_format("log10", scales::math_format(10^.x))
 ) +
  annotation_logticks(side = "b", outside=TRUE) +
  coord_cartesian(clip = "off") + 
  xlab("Number of bugs") + 
  ylab("Density")

Visually we seem to have a very good fit.

2.3.3 Diagnostics

Divergent transitions:

np <- nuts_params(m2)
sum(subset(np, Parameter == "divergent__")$Value)

## [1] 0

\(\widehat{R}\) ratio:

max(rhat(m2), na.rm = TRUE)

## [1] 1.004551

Effective sample size (ESS):

min(neff_ratio(m2), na.rm = TRUE)

## [1] 0.1137986

Finally, we have also visually inspected the trace plots, but since this model estimates more than \(1787\) parameters we choose not to show it here.

In the end, the model passed all diagnostics, but before we continue analyzing the output of this, more complex model, we might as well see which of our three models, \(\{\mathcal{M}_0, \mathcal{M}_1, \mathcal{M}_2\}\), has the best out of sample prediction capabilities.

3 Model comparison

In Bayesian statistics it has been quite common to compare models’ relative out of sample prediction capabilities. This has often been done by using information theoretical concepts like, e.g., AIC, BIC, DIC, and WAIC. Up until recently, many would argue that WAIC (Widely Applicable Information Criterion, or Watanabe-Akaike Information Criterion) was state of the art concerning model comparison. However, by standing on the shoulder of giants, and all that, things have changed.

Using Pareto-smoothed importance sampling (PSIS), together with leave-one-out cross validation (LOO), we can do all the things that one could do with WAIC, but with the added benefit of having diagnostics available that tells us when things go wrong (Vehtari, Gelman, and Gabry 2017).

First, we add the information criteria, i.e., LOO, to the objects.

m0 <- add_criterion(m0, criterion = "loo")
m1 <- add_criterion(m1, criterion = "loo")
m2 <- add_criterion(m2, criterion = "loo")

Next, we calculate the relative out of sample prediction capabilities of the models.

(lc <- loo(m0,m1,m2, moment_match = TRUE))

## Output of model 'm0':
## 
## Computed from 10000 by 1127 log-likelihood matrix
## 
##          Estimate    SE
## elpd_loo  -7289.8  76.8
## p_loo        52.7  18.5
## looic     14579.6 153.7
## ------
## Monte Carlo SE of elpd_loo is NA.
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. n_eff
## (-Inf, 0.5]   (good)     1121  99.5%   920       
##  (0.5, 0.7]   (ok)          2   0.2%   219       
##    (0.7, 1]   (bad)         2   0.2%   45        
##    (1, Inf)   (very bad)    2   0.2%   3         
## See help('pareto-k-diagnostic') for details.
## 
## Output of model 'm1':
## 
## Computed from 4000 by 1127 log-likelihood matrix
## 
##          Estimate    SE
## elpd_loo  -5496.0  60.2
## p_loo        22.6   2.2
## looic     10992.0 120.4
## ------
## Monte Carlo SE of elpd_loo is 0.1.
## 
## All Pareto k estimates are good (k < 0.5).
## See help('pareto-k-diagnostic') for details.
## 
## Output of model 'm2':
## 
## Computed from 4000 by 1127 log-likelihood matrix
## 
##          Estimate    SE
## elpd_loo  -5325.7  59.5
## p_loo       397.5  15.7
## looic     10651.4 118.9
## ------
## Monte Carlo SE of elpd_loo is NA.
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. n_eff
## (-Inf, 0.5]   (good)     652   57.9%   584       
##  (0.5, 0.7]   (ok)       355   31.5%   77        
##    (0.7, 1]   (bad)      110    9.8%   30        
##    (1, Inf)   (very bad)  10    0.9%   5         
## See help('pareto-k-diagnostic') for details.
## 
## Model comparisons:
##    elpd_diff se_diff
## m2     0.0       0.0
## m1  -170.3      23.8
## m0 -1964.1      58.4

We see that there’s some diagnostics showing that we have Pareto \(k>0.7\) in \(\mathcal{M}_0\) and \(\mathcal{M}_2\). In this particular case, running \(10\)-fold cross validation, one will see that the estimations of the LOO information criteria is stable, i.e., the criteria is precisely estimated using moment matching correction to the importance sampling for the problematic observations (Paananen et al. 2019).

At the bottom of the output, we see that \(\mathcal{M}_2\) seems to take the lead (placed on the first row). If we assume \(z_{99\%} = 2.58\), then one clearly sees that it’s considerably better, i.e., \(\mathrm{CI}_{95\%}\)[-231.7, -108.9] does not cross zero (after all it’s more than \(7\) SE away from \(\mathcal{M}_1\)). In short, \(\mathcal{M}_2\) is very likely the best model we have in our arsenal at the moment, hence, we declare it to be \(\mathcal{M}\),

M <- m2

4 Inferences on \(\mathcal{M}\)

A violin plot summarizes the fit for each language (\(y\) is our empirical data, and \(y_{\mathrm{rep}}\) are draws from the posterior probability distribution):

pp_check(M, type = "violin_grouped", group = "language_id", y_draw="points") +
  scale_y_continuous(trans = "log2") +
  scale_x_discrete(labels=levels(fse.data$language)) +
  theme(axis.text.x = element_text(angle=45, hjust=1))

We can also check the conditional effects on our population-level parameters.

For the first three effects we clearly see something is happening when moving along the \(x\)-axis. With the last effect, ‘insertions,’ it is not as clear; mainly due to the 95% credible interval that dramatically increases. The latter indicates that our variable ‘insertions’ might not be suitable for inclusion in the model.

What are our population-level estimates? We refrain from printing the group effects (random effects) since it will take up a lot of space but they can be found in the appendix.⁸

fixef(M)

##                       Estimate  Est.Error        Q2.5       Q97.5
## Intercept          -2.02125595 0.17071437 -2.35306603 -1.67801247
## devs_log            0.06434463 0.01841558  0.02835247  0.10157985
## max_commit_age_log  0.06317768 0.02036218  0.02170501  0.10266619
## commits_log         0.96848155 0.02288104  0.91905692  1.01092981
## insertions_log      0.02440615 0.01819853 -0.01046735  0.06153201

In short, three of the \(\beta\) estimates are significant on 95%, but ‘insertions’ is not.

Let’s next plot our group effects for ‘language’ (varying according to our slopes ‘devs,’ ‘max_commit_age,’ ‘commits,’ and ‘insertions’). The inner interval is set to 50% of the probability mass (light blue), while the outer interval is set to 90%.

Generally speaking, the slope ‘insertions’ shows more variance, as does, but slightly harder to see, ‘commits’ (they both have \(\sigma=0.05\)). It’s not visible in the plots (but when we look at the summary()), the strongest correlation is that ‘commits’ is negatively correlated with ‘insertions’ (\(r=0.46\)). Once again we see ‘insertions’ popping up, and in the above plot Typescript looks suspicious, i.e., something seems to be a bit fishy.

4.1 Projects’ variability

By plotting our slopes for, e.g.,‘insertions,’ we could see some suspect things. The variability is larger compared to other slopes, but also that Typescript sticks out. This we all got by designing a model that implemented varying intercepts (according to language) with varying slopes (according to the population-level effects insertions, commits, max commit age, and number of developers). However, we potentially have another variable we could use as a varying intercept, i.e., the project id.

Recall that in \(\mathcal{M}_1\) we did not have project id as a varying intercept, and comparing that model to our final model, we saw that \(\mathcal{M}_1\) simply could not compete in out of sample prediction. In short, adding project id clearly made a large difference from an information theoretical point of view.

Examining project id we see that there are 729 projects. Let’s see the distribution of rows for each project, i.e., how many rows do we have for each project?

A clear majority of the projects only have one row. How peculiar, dynamic Hamiltonian Monte Carlo (HMC) still cuts through this like a warm knife through butter. Let’s plot the prior probability distribution, \(\mathrm{Weibull}(2,1)\), and compare it with our posterior.

There’s no question that the data has swamped our prior. One clearly sees our wide prior (dark) compared to the ‘spike,’ which is our posterior.

In fact, what we’re seeing is the effect of multi-level modeling taking place. In order for us to avoid overfitting, projects with more data will inform projects with little data. It’s a remarkable fact that HMC can sample this so nicely in the end. Making use of the project id allowed us to improve the out of sample predictions considerably!

Let’s plot a random sample of \(10\) projects to get a feeling for the variability among projects (this plot will of course look different every time we plot it). We’ll plot the \(x\)-axis using logarithm base \(10\).

So, projects contain a considerable within and between variability, which is most likely why the model makes such a jump in out of sample predictions when we add projects as a varying intercept. More than this we don’t need to know, after all, what is a significant project? We’d rather compare the variability between and withing languages.

4.2 Predictions

As a first step, we can now set our predictors to the original median values and see how the languages differ. Let’s query our posterior probability distribution about the median bugs for each language (all things being equal).

n_bugs	lang
79	Objective-C
78	C++
74	C
73	Php
72	Coffeescript
70	Python
68	Javascript
66	Erlang
65	C#
65	Java
63	Ruby
61	Perl
60	Go
55	Haskell
54	Scala
52	Clojure
51	Typescript

The table only provides us with a rough overview of language differences when setting continuous covariates to their medians. We should try using our covariates at different levels and plot the outcome with uncertainty to see how it varies depending on the level.

The plots below shows posterior predictions of the number of bugs, when our covariates (predictors) are set to their empirical max, median, and minimum levels (note that the \(y\)-axes use logarithm base \(10\)!) The max level clearly shows uncertainty, but that uncertainty is even more pronounced for the other two levels.

What is clearly evident is that the order (languages are plotted in a descending order from left to right, where the right-most languages have fewer bugs) changes depending on covariates’ settings. So there’s no black or white here, i.e., depending on the covariates’ values different languages perform differently.

Which of the above plots tells the truth? Well, the truth is connected to your reality, i.e., practical significance. Hence, it varies a lot depending on contextual factors. In our case these factors are represented by the covariates’ levels/settings, but many other such factors exist in the form of omitted variable bias.

4.2.1 Setting our own reality

Let’s assume that we have a project with \(30\) developers which will run for two years using Ruby and Python. Over these two years we expect on average \(2\) commits per day (over 365 days) \(2*365\), max commit age of \(30*2*365\), and insertions ending up around \(30*10*2*365\)

newdata <- data.frame(language_id=c(14,15), # Python and Ruby
                      devs_log=log(30),
                      project_id = NA,
                      max_commit_age_log=log(2*365),
                      commits_log=log(30*2*365),
                      insertions_log=log(30*10*2*365))

PPC_new <- posterior_predict(M, newdata=newdata)
names <- c("Python", "Ruby")
colnames(PPC_new) <- names
PPC_new <- as.data.frame(PPC_new) %>%
                  gather(lang, bugs)

ggplot(PPC_new, aes(x=lang, y=bugs)) +
  geom_violin() +
  xlab("") +
  ggtitle("Custom (developers = 30, max commit age = 2*365, \ncommits = 30*2*365, insertions = 30*10*2*365)")

All things being equal, you should consider using Ruby in your project. Ruby has a shorter tail (i.e., has a maximum number of bugs lower than Python). Additionally, the bulk of the probability mass is lower than for Python.

But do we have a ‘statistically significant difference’ between these two languages, or for that matter among any of the languages?

4.3 Effect sizes

We’ll use covariate values from the empirical data, \(y\), as input, with the idea that the sample is representative of the population. This is really not needed, i.e., we can set the values to basically whatever we want, but this way it’ll be easier to conduct inferences for all combinations of languages, and since we have a posterior probability distribution we usually don’t have to worry about correcting for multiple comparisons.

If we use a credible level of 0.999 (i.e., 0.1%) we have 116 pairwise differences that are of interest. Then, for each contrast we calculate the 95% Highest Posterior Density Interval (HPDI).⁹ If the HPDI is strictly positive or negative (i.e., does not cover 0) then we have a significant difference on the 95% HPDI concerning the difference of two languages.

First, we calculate the contrast between each language (136 combinations when we have 17 languages). Then it’s very easy to check the distribution of the difference between any two languages and in that way also investigate effect sizes (there’s no need to look at point estimates since we have a probability distribution of the differences between two languages, with 95% highest posterior density interval).

If we plot the difference between C# and C (above) we see that it’s positive, indicating that C almost always performs worse than C# (if it would have been negative then C# would’ve been worse).

Let’s also look (below) at cases where there is no significant difference (Coffeescript vs. Go)

We see that it crosses 0 on the \(x\)-axis, and sometimes Go actually is worse on average.

How often is sometimes?

es <- data.frame(table(sign(contrasts[[60]])))
(foo <- 1 - (es[1,2] /(es[1,2] + es[2,2])))

## [1] 0.083

Well, in this case 8.3% of the time.

Finally, here’s an example when there is a clear difference in the other direction (Objective-C vs. Ruby), i.e., when the values are clearly negative and, so, Objective-C is clearly worse.

Once again, these are the distributions in differences that we get when we reuse the covariate statistics from our sample and average over them. As we saw previously, with other covariates we get other outcomes. Rarely does a statistical analysis provide us with a dichotomous answer, i.e., the ‘truth’ is seldom that simple.

5 Bayesian variable selection for \(\textrm{NB}(\lambda,\phi)\)

We did not include this step in the original manuscript since one of the purposes with the study was to also use the same independent variables (i.e., predictors) as in the previously published work. However, we would argue that this step should be used if one wants to gain robust out of sample predictions.

Here we’ll conduct variable selection for negative-binomial (NB). Unfortunately, the Stan team does not yet (January 2021) support variable selection for the NB.¹⁰ However, there’s a package that uses MCMC for variable selection of NB distributions by Dvorzak and Wagner (2016).

Set our outcome variable:

y <- fse.data$n_bugs

Select response variables and store them as a matrix. We remove ‘n_bugs’ (our \(y\)), original variables not transformed, and the factor variables ‘language’ and ‘project’ (since they are indicators we’ll use later as varying intercepts).

X <- as.matrix(fse.data[,c(-1:-8, -13)])

Run the variable selection,

varsel_res <- negbinBvs(y = y, X = X)

and plot the chain and output the summary:

plot(varsel_res, burnin = FALSE, thin = FALSE)

summary(varsel_res)

## Bayesian variable selection for the negative binomial model:
## 
## Call:
## negbinBvs(y = y, X = X)
## 
## 
## MCMC:
## M = 8000 draws after a burn-in of 2000
## BVS started after 1000 iterations
## Thinning parameter: 1
## 
## Acceptance rate for rho:
## 69.46%
## 
## Prior: spike-and-slab prior with Student-t slab [V=5]
## 
## b0[0] b0[1] b0[2] b0[3] b0[4] 
##     0     0     0     0     0 
## w[a] w[b]   c0   C0 
##    1    1    2    1 
## 
## 
## Model averaged posterior means, estimated posterior inclusion
##  probabilities and 95%-HPD intervals:
## 
##             Estimate P(.=1) 95%-HPD[l] 95%-HPD[u]
## (Intercept)   -1.954     NA     -2.165     -1.760
## beta.1         0.982  1.000      0.956      1.010
## beta.2         0.016  0.649      0.000      0.038
## beta.3         0.061  1.000      0.029      0.092
## beta.4         0.072  1.000      0.043      0.098
## rho            6.340     NA      5.741      6.973

\(\beta_2\) \(\mathrm{P}(.= 1)\) is \(\approx 0.5\) so we should think carefully about including it, i.e., the variable ‘insertions.’ We also have fairly strong argument for ‘max_commit_age’ being a bit shaky (\(\beta_3\) in the above case) by looking at the lower 95% highest posterior density (95%-HPD[l]). For the sake of completeness, and to follow the previous studies, we included them anyway (in our analysis above). But, if we would follow a principled Bayesian workflow, then many would argue that we should not use the predictor ‘insertions.’ It’s also clear, by looking at our analysis further above, that ‘insertions’ is not significant and has questionable contributions to the prediction capabilities of the model (i.e., it contributes with uncertainty, and not much else).

A model without the variable ‘insertions’ would look like this,

m_no_insertions <- brm(n_bugs ~ 1 + devs_log + max_commit_age_log + commits_log +
                         (1 + devs_log + max_commit_age_log + 
                            commits_log | language_id) +
                         (1 | project_id),
                       family = negbinomial,
                       data = fse.data,
                       prior = p,
                       control = list(adapt_delta = 0.95)
                       )

In the model above we’ve set weakly informative priors and the model’s diagnostics indicate that we’ve converged towards a stationary posterior.

However, if we compare the parameter estimates from this model to \(\mathcal{M}\), we see that there are very little changes.

fixef(M)

##                       Estimate  Est.Error        Q2.5       Q97.5
## Intercept          -2.02125595 0.17071437 -2.35306603 -1.67801247
## devs_log            0.06434463 0.01841558  0.02835247  0.10157985
## max_commit_age_log  0.06317768 0.02036218  0.02170501  0.10266619
## commits_log         0.96848155 0.02288104  0.91905692  1.01092981
## insertions_log      0.02440615 0.01819853 -0.01046735  0.06153201

fixef(m_no_insertions)

##                       Estimate  Est.Error        Q2.5       Q97.5
## Intercept          -1.90272951 0.13227542 -2.16282079 -1.63849895
## devs_log            0.06199018 0.01901210  0.02405677  0.10034399
## max_commit_age_log  0.05612958 0.02109243  0.01487555  0.09717222
## commits_log         0.99591989 0.01499997  0.96447794  1.02326755

Let’s now turn our attention to the TOPLAS dataset, which is a ‘cleaned’ subset of the FSE dataset. The above analysis will now guide us when we analyze the TOPLAS dataset. We have followed the same workflow as above for the TOPLAS analysis, but we do not report on all steps.

6 Computational environment

sessionInfo()

## R version 4.0.3 (2020-10-10)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur 10.16
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## attached base packages:
## [1] parallel  stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] ggridges_0.5.3          patchwork_1.1.1         latex2exp_0.4.0        
##  [4] LaplacesDemon_16.1.4    bayesplot_1.8.0         ggthemes_4.2.4         
##  [7] rethinking_2.13         rstan_2.26.0.9000       StanHeaders_2.26.0.9000
## [10] brms_2.14.6             Rcpp_1.0.6              pogit_1.2.0            
## [13] ggplot2_3.3.3           kableExtra_1.3.1        tidyr_1.1.2            
## [16] dplyr_1.0.3            
## 
## loaded via a namespace (and not attached):
##   [1] minqa_1.2.4          colorspace_2.0-0     ellipsis_0.3.1      
##   [4] rsconnect_0.8.16     markdown_1.1         base64enc_0.1-3     
##   [7] rstudioapi_0.13      farver_2.0.3         DT_0.17             
##  [10] fansi_0.4.2          mvtnorm_1.1-1        xml2_1.3.2          
##  [13] codetools_0.2-18     bridgesampling_1.0-0 splines_4.0.3       
##  [16] knitr_1.31           shinythemes_1.2.0    projpred_2.0.2      
##  [19] jsonlite_1.7.2       nloptr_1.2.2.2       shiny_1.6.0         
##  [22] compiler_4.0.3       httr_1.4.2           backports_1.2.1     
##  [25] assertthat_0.2.1     Matrix_1.3-2         fastmap_1.1.0       
##  [28] cli_2.2.0            later_1.1.0.1        htmltools_0.5.1.1   
##  [31] prettyunits_1.1.1    tools_4.0.3          igraph_1.2.6        
##  [34] coda_0.19-4          gtable_0.3.0         glue_1.4.2          
##  [37] reshape2_1.4.4       V8_3.4.0             vctrs_0.3.6         
##  [40] nlme_3.1-151         crosstalk_1.1.1      xfun_0.20           
##  [43] stringr_1.4.0        ps_1.5.0             lme4_1.1-26         
##  [46] rvest_0.3.6          mime_0.9             miniUI_0.1.1.1      
##  [49] lifecycle_0.2.0      gtools_3.8.2         statmod_1.4.35      
##  [52] MASS_7.3-53          zoo_1.8-8            scales_1.1.1        
##  [55] colourpicker_1.1.0   promises_1.1.1       Brobdingnag_1.2-6   
##  [58] inline_0.3.17        shinystan_2.5.0      gamm4_0.2-6         
##  [61] yaml_2.2.1           curl_4.3             gridExtra_2.3       
##  [64] loo_2.4.1            stringi_1.5.3        highr_0.8           
##  [67] dygraphs_1.1.1.6     checkmate_2.0.0      boot_1.3-26         
##  [70] pkgbuild_1.2.0       shape_1.4.5          rlang_0.4.10        
##  [73] pkgconfig_2.0.3      matrixStats_0.57.0   evaluate_0.14       
##  [76] lattice_0.20-41      purrr_0.3.4          labeling_0.4.2      
##  [79] rstantools_2.1.1     htmlwidgets_1.5.3    tidyselect_1.1.0    
##  [82] processx_3.4.5       plyr_1.8.6           magrittr_2.0.1      
##  [85] bookdown_0.21        R6_2.5.0             generics_0.1.0      
##  [88] DBI_1.1.1            pillar_1.4.7         withr_2.4.1         
##  [91] mgcv_1.8-33          xts_0.12.1           abind_1.4-5         
##  [94] tibble_3.0.5         crayon_1.3.4         utf8_1.1.4          
##  [97] rmarkdown_2.6        grid_4.0.3           callr_3.5.1         
## [100] threejs_0.3.3        digest_0.6.27        webshot_0.5.2       
## [103] xtable_1.8-4         httpuv_1.5.5         RcppParallel_5.0.2  
## [106] stats4_4.0.3         munsell_0.5.0        viridisLite_0.3.0   
## [109] shinyjs_2.0.0

end.time <- Sys.time()
round((end.time - start.time), 3)

## Time difference of 21.931 mins

7 Appendix

Estimated random effects from the candidate model \(\mathcal{M}\). The languages (\(1\)–\(17\)) are the following,

levels(fse.data$language)

##  [1] "C"            "C#"           "C++"          "Clojure"      "Coffeescript"
##  [6] "Erlang"       "Go"           "Haskell"      "Java"         "Javascript"  
## [11] "Objective-C"  "Perl"         "Php"          "Python"       "Ruby"        
## [16] "Scala"        "Typescript"

ranef(M)

## $language_id
## , , Intercept
## 
##       Estimate Est.Error       Q2.5     Q97.5
## 1   0.16732217 0.2520515 -0.3416274 0.6553464
## 2  -0.25526048 0.2927475 -0.8666162 0.2939253
## 3   0.12293541 0.2748935 -0.3930773 0.6776701
## 4  -0.27421410 0.3509270 -0.9901859 0.3824851
## 5   0.14428476 0.2891762 -0.4240638 0.7357660
## 6   0.01706006 0.3409835 -0.6887694 0.6698310
## 7  -0.04542180 0.3233300 -0.6955682 0.5838824
## 8  -0.02241414 0.3631452 -0.7462205 0.6824683
## 9  -0.03927441 0.2890932 -0.5936530 0.5440560
## 10 -0.02058118 0.2132276 -0.4396024 0.3944612
## 11 -0.36018404 0.3355958 -1.0833555 0.2529063
## 12 -0.38296661 0.3465104 -1.0690746 0.2960590
## 13  0.07705889 0.2795713 -0.4633927 0.6334283
## 14 -0.12269316 0.2432051 -0.6229433 0.3385126
## 15 -0.12970314 0.2749988 -0.6999160 0.4084258
## 16 -0.26972829 0.3594526 -1.0020859 0.4371881
## 17  1.34270611 0.3894801  0.5880608 2.1444870
## 
## , , devs_log
## 
##         Estimate  Est.Error        Q2.5      Q97.5
## 1   0.0189667130 0.03037167 -0.03480087 0.08624825
## 2   0.0244181708 0.03195438 -0.02797041 0.09849048
## 3  -0.0381213537 0.03199088 -0.10983537 0.01302325
## 4   0.0042794428 0.03202971 -0.06118176 0.07237525
## 5  -0.0164103955 0.03125846 -0.08880505 0.03738975
## 6   0.0202454753 0.03786457 -0.04574336 0.10729878
## 7   0.0090776847 0.03073040 -0.05360626 0.07517994
## 8   0.0078297371 0.03005337 -0.05267532 0.07310074
## 9   0.0151285514 0.02915377 -0.03699645 0.07963961
## 10 -0.0046694195 0.02342574 -0.05273010 0.04025200
## 11  0.0121475984 0.03293982 -0.05095713 0.08323022
## 12 -0.0003249481 0.03393086 -0.07406130 0.06672694
## 13 -0.0059970647 0.02981003 -0.07055650 0.05328807
## 14 -0.0100390420 0.02758922 -0.07195053 0.03978134
## 15 -0.0068515561 0.02830826 -0.06824258 0.04912273
## 16  0.0030362798 0.03152732 -0.06330610 0.06658342
## 17 -0.0344769394 0.04230745 -0.12598273 0.03954601
## 
## , , max_commit_age_log
## 
##         Estimate  Est.Error        Q2.5      Q97.5
## 1   0.0180530857 0.02888718 -0.03080884 0.08498319
## 2   0.0187590458 0.03158799 -0.03310582 0.09539045
## 3  -0.0065507082 0.02818098 -0.06617946 0.04656438
## 4  -0.0062533394 0.03120212 -0.06978971 0.05811904
## 5  -0.0017805693 0.02804112 -0.05969508 0.05985541
## 6   0.0204820036 0.03605555 -0.03699768 0.10766476
## 7   0.0015023690 0.02940894 -0.05931314 0.06243409
## 8  -0.0019115318 0.03067204 -0.06625945 0.05896649
## 9  -0.0103107779 0.02855185 -0.07453784 0.04029509
## 10  0.0005265331 0.02249269 -0.04480085 0.04712385
## 11  0.0185836985 0.03319438 -0.03797210 0.09583997
## 12 -0.0098933882 0.03282367 -0.08618724 0.05006150
## 13  0.0044603918 0.02830672 -0.05451377 0.06586291
## 14  0.0113129136 0.02659867 -0.03720003 0.07326270
## 15  0.0036112177 0.02849275 -0.05285337 0.06625733
## 16 -0.0191509244 0.03739848 -0.11408039 0.04051266
## 17 -0.0384225426 0.04297382 -0.13663183 0.03096799
## 
## , , commits_log
## 
##        Estimate  Est.Error         Q2.5      Q97.5
## 1   0.008337577 0.02907992 -0.048140827 0.06972154
## 2  -0.033178004 0.04020730 -0.124063007 0.03274397
## 3   0.018931405 0.03387065 -0.042923461 0.09154683
## 4   0.005559826 0.04562738 -0.083590098 0.10226049
## 5  -0.003486763 0.04284132 -0.093054337 0.08725170
## 6  -0.039558591 0.04734672 -0.145823949 0.03876512
## 7  -0.005245115 0.04200152 -0.091274688 0.07986839
## 8  -0.018214258 0.04561665 -0.124412665 0.06282227
## 9   0.025328048 0.03980472 -0.045303781 0.11386641
## 10  0.045185130 0.03306233 -0.008183022 0.11815396
## 11 -0.010897630 0.04242917 -0.099749409 0.07336005
## 12 -0.005616153 0.04473795 -0.098318876 0.08688512
## 13 -0.003223747 0.03910627 -0.084329045 0.07550462
## 14  0.014730311 0.03235058 -0.043917648 0.08713381
## 15 -0.042098969 0.04268690 -0.137951217 0.02684256
## 16  0.003495511 0.04380309 -0.084590231 0.09608489
## 17  0.054980884 0.05255546 -0.037263902 0.16596516
## 
## , , insertions_log
## 
##        Estimate  Est.Error        Q2.5       Q97.5
## 1  -0.023266020 0.02449626 -0.07313556  0.02238382
## 2   0.024002078 0.03285616 -0.03738884  0.09326814
## 3   0.011625255 0.02742197 -0.04230698  0.06724491
## 4   0.005613337 0.04095666 -0.07424843  0.08721006
## 5   0.004184757 0.03717211 -0.06772229  0.07979320
## 6   0.002753943 0.03679679 -0.06691301  0.07720890
## 7  -0.003619793 0.03850586 -0.07811510  0.07692604
## 8  -0.003866558 0.04401515 -0.08872699  0.08581301
## 9  -0.006394690 0.03327855 -0.07742165  0.05425693
## 10 -0.016414272 0.02377136 -0.06560131  0.02933206
## 11  0.044871998 0.03726016 -0.02681116  0.12499879
## 12  0.041261094 0.04014471 -0.03392014  0.12497978
## 13  0.004638859 0.03227035 -0.06035085  0.07231875
## 14  0.007353516 0.02653848 -0.04643870  0.05820578
## 15  0.032947559 0.03607074 -0.03301834  0.11031232
## 16  0.019421399 0.04011481 -0.05793193  0.10336349
## 17 -0.144898680 0.03918707 -0.22352252 -0.06980498
## 
## 
## $project_id
## , , Intercept
## 
##          Estimate  Est.Error         Q2.5         Q97.5
## 1   -0.7868424120 0.24537843 -1.274469048 -0.3074238617
## 2    0.1237329896 0.15469684 -0.177480358  0.4252849197
## 3   -0.1581682384 0.18539854 -0.526837390  0.2069117005
## 4   -0.1608717734 0.21894865 -0.597075361  0.2791943034
## 5    0.2091952793 0.18251528 -0.142444687  0.5677284942
## 6   -0.2107026345 0.18771880 -0.572773442  0.1648760713
## 7   -0.0860571908 0.14674596 -0.371474253  0.2013441209
## 8    0.2410415237 0.21003998 -0.162500861  0.6472721571
## 9   -0.0087570243 0.18754639 -0.375223480  0.3543228099
## 10   0.1706063310 0.18177459 -0.184053653  0.5353254422
## 11   0.0851548155 0.17409174 -0.253800061  0.4382993814
## 12  -0.0841914895 0.18808500 -0.449244344  0.2884992013
## 13  -0.1964179035 0.13455680 -0.458263727  0.0652422160
## 14  -0.0577547919 0.17974933 -0.410135095  0.3017591758
## 15  -0.1827873065 0.19383901 -0.579594369  0.2027938110
## 16  -0.0078516487 0.16256009 -0.335184911  0.3108204515
## 17   0.0022061172 0.23764936 -0.474314736  0.4698097135
## 18   0.0753981298 0.14704702 -0.207473908  0.3687161411
## 19  -0.2131343755 0.18549329 -0.571323412  0.1581826382
## 20  -0.2288066024 0.16535812 -0.546946004  0.0979644533
## 21   0.0810660119 0.18605361 -0.281242971  0.4496716991
## 22   0.1533940905 0.20769073 -0.249272238  0.5579066317
## 23   0.2031738895 0.16965600 -0.127286595  0.5441149639
## 24  -0.0546316066 0.19683602 -0.442028757  0.3203285366
## 25   0.0262418681 0.17385609 -0.311223115  0.3707863935
## 26  -0.3749923252 0.20071474 -0.759333322  0.0132812330
## 27   0.1843765451 0.17376474 -0.162833270  0.5195005943
## 28  -0.1251702962 0.15166983 -0.424313620  0.1742838549
## 29  -0.0127376422 0.15972369 -0.319106258  0.3019298612
## 30  -0.0724613396 0.16985570 -0.407072844  0.2653522406
## 31  -0.1618326740 0.21571570 -0.583491952  0.2612018861
## 32   0.1177548897 0.21008974 -0.300233159  0.5410338468
## 33   0.0886886014 0.15406587 -0.206903600  0.3900428986
## 34   0.0511474143 0.12786863 -0.195578746  0.3067866856
## 35   0.0874307741 0.17077627 -0.232627892  0.4167381792
## 36  -0.1061835252 0.18188138 -0.452633777  0.2523986448
## 37  -0.1259277335 0.20018718 -0.509252905  0.2662135543
## 38   0.3089614688 0.13298915  0.046167894  0.5602947065
## 39   0.7084940116 0.14977796  0.417646209  0.9950659407
## 40   0.0090069818 0.17463163 -0.331054654  0.3449804970
## 41   0.0918141757 0.17858995 -0.258014143  0.4400172001
## 42  -0.1970481910 0.19270621 -0.570188806  0.1830124327
## 43  -0.0258489639 0.17227331 -0.349901496  0.3211736230
## 44   0.0354768413 0.19247594 -0.345824504  0.4121763874
## 45  -0.1211607664 0.23357154 -0.586096123  0.3204526378
## 46  -0.0597882541 0.14106567 -0.336328344  0.2243957574
## 47  -0.0994905267 0.18326449 -0.453348008  0.2679187709
## 48   0.2619270136 0.17046891 -0.071927398  0.6041714908
## 49   0.0196975818 0.18680286 -0.343937949  0.3920128130
## 50  -0.0203777532 0.14918655 -0.308864770  0.2694322559
## 51  -0.0557885642 0.14917520 -0.342577866  0.2518351541
## 52  -0.0129911587 0.15970277 -0.321005373  0.3048783357
## 53  -0.0559535402 0.19612973 -0.437474612  0.3152743481
## 54  -0.0664309352 0.17197447 -0.401365684  0.2782614195
## 55  -0.3156562758 0.13041910 -0.561780691 -0.0491327818
## 56  -0.0851509341 0.17737950 -0.424562095  0.2623267481
## 57  -0.0254479486 0.23023337 -0.475001699  0.4169644549
## 58  -0.0498373908 0.17318650 -0.393734392  0.2976832678
## 59  -0.0939079380 0.16378946 -0.415012066  0.2364687495
## 60   0.0553925839 0.17110953 -0.271135506  0.3957462050
## 61   0.1527314302 0.16638052 -0.172213140  0.4826141890
## 62  -0.2723981682 0.21897170 -0.700071349  0.1424541356
## 63   0.1777944417 0.17265694 -0.149349365  0.5304222380
## 64  -0.0836936092 0.17536144 -0.421486124  0.2625501458
## 65  -0.0627006131 0.18254394 -0.418165061  0.2940935311
## 66  -0.3024841767 0.18674173 -0.654887687  0.0777130567
## 67   0.0325770216 0.22575792 -0.412737262  0.4691490476
## 68   0.7058621018 0.16871580  0.374062718  1.0355961676
## 69  -0.1667732571 0.15051992 -0.458388074  0.1316181349
## 70  -0.1292679183 0.22190854 -0.565363591  0.3101723007
## 71  -0.1575354989 0.19692172 -0.541111124  0.2344538793
## 72   0.0892771814 0.17350994 -0.251178911  0.4279755381
## 73   0.1629781521 0.17607728 -0.178236957  0.5153008722
## 74  -0.0080217959 0.17916891 -0.357007742  0.3421294543
## 75   0.1223547988 0.19659732 -0.263683673  0.5120137058
## 76  -0.0341472933 0.21943346 -0.477402546  0.3847388468
## 77  -0.1459731442 0.17866611 -0.497310292  0.1942532637
## 78  -0.0909423767 0.17727704 -0.430500495  0.2613014382
## 79  -0.1169540590 0.15497117 -0.418096816  0.1972618070
## 80  -0.1801177867 0.22712347 -0.623411843  0.2560920329
## 81  -0.4619003799 0.18698540 -0.825717781 -0.0900517481
## 82  -0.0808186175 0.18169217 -0.425126046  0.2730998150
## 83   0.0997700945 0.20218681 -0.292023824  0.5097445697
## 84  -0.2765361818 0.20934613 -0.688383452  0.1308736201
## 85   0.0689258917 0.17778350 -0.271116577  0.4168764536
## 86  -0.1518061167 0.14968458 -0.444133475  0.1521801541
## 87  -0.0675881867 0.22635931 -0.512601921  0.3717493572
## 88  -0.1503216227 0.14652785 -0.436518290  0.1349140372
## 89   0.1340894642 0.20195480 -0.261227867  0.5249870516
## 90  -0.1512947457 0.21487735 -0.569588717  0.2619033921
## 91  -0.0979044185 0.18201985 -0.463432168  0.2580192408
## 92  -0.1688075909 0.13312637 -0.424243557  0.0928254479
## 93  -0.1515296711 0.20170848 -0.540819353  0.2482180758
## 94  -0.0259236391 0.17469490 -0.367894313  0.3152722380
## 95  -0.4990225601 0.18440621 -0.865963825 -0.1398937481
## 96   0.0178517964 0.17363395 -0.323024125  0.3586340209
## 97  -0.0085278734 0.18009480 -0.353845421  0.3606363697
## 98  -0.0325340138 0.22713688 -0.478971832  0.4132309419
## 99   0.1379324460 0.14803592 -0.145191505  0.4332853561
## 100  0.1653339941 0.17553168 -0.177031871  0.5219679786
## 101  0.1364136978 0.18306353 -0.215592544  0.5042676252
## 102 -0.0397899137 0.15552702 -0.338511975  0.2737284782
## 103  0.4537883314 0.09279979  0.273419368  0.6406365920
## 104  0.0912783379 0.19529116 -0.282302651  0.4898144974
## 105 -0.1795509392 0.15702903 -0.486838489  0.1390593444
## 106 -0.0513925225 0.19041959 -0.416756452  0.3269138236
## 107  0.0474123890 0.12012084 -0.186012663  0.2854042529
## 108  0.0279176547 0.17442068 -0.312082168  0.3755652555
## 109  0.2527094983 0.14286465 -0.019346950  0.5417861707
## 110  0.0421470644 0.17145342 -0.284140050  0.3929283896
## 111 -0.4241297815 0.17210477 -0.753118105 -0.0778067075
## 112 -0.1391589812 0.18376125 -0.489418578  0.2186137510
## 113  0.1105344681 0.17042828 -0.221095431  0.4590505012
## 114 -0.2773043036 0.19233518 -0.649837866  0.1117379254
## 115 -0.1019711557 0.22866606 -0.558956263  0.3442901469
## 116 -0.0704007442 0.13952727 -0.337848526  0.2077344320
## 117  0.2080644318 0.18489996 -0.156096220  0.5783376266
## 118  0.0383882225 0.18167277 -0.322056571  0.3976871676
## 119  0.2023900341 0.11486019 -0.024256351  0.4250822976
## 120  0.0136117268 0.16423867 -0.299516902  0.3486093797
## 121  0.0848807026 0.17623908 -0.258163488  0.4383660872
## 122  0.1123803119 0.14214747 -0.164872344  0.3959430387
## 123 -0.0306344959 0.18329581 -0.386868590  0.3332194676
## 124 -0.2468873380 0.17075566 -0.572459103  0.0867521382
## 125  0.3370556635 0.16142630  0.026155499  0.6572349812
## 126 -0.1652395058 0.23568849 -0.636448573  0.3027992667
## 127 -0.0511445008 0.21808175 -0.482561348  0.3804931436
## 128  0.1670815050 0.18983057 -0.206617794  0.5429997812
## 129 -0.0234347818 0.17920604 -0.372555866  0.3373904849
## 130  0.1052508412 0.17994449 -0.242055786  0.4648717204
## 131 -0.0812402089 0.17790888 -0.424837225  0.2775649831
## 132 -0.0602468803 0.20219678 -0.458717225  0.3381778868
## 133 -0.0597571524 0.15681363 -0.360800864  0.2474338210
## 134 -0.0685710057 0.16180222 -0.379873611  0.2510393293
## 135  0.0468359169 0.21780016 -0.383496371  0.4760215357
## 136  0.1307348121 0.19088394 -0.244290847  0.4981839722
## 137  0.3258250054 0.13505680  0.073620645  0.5919872443
## 138 -0.4466488399 0.21173474 -0.853113325 -0.0350436580
## 139 -0.3499781878 0.19376021 -0.721789980  0.0338187433
## 140 -0.1405271914 0.17751642 -0.481528368  0.2111777917
## 141 -0.1023784975 0.14415837 -0.377907196  0.1877461623
## 142  0.0923617435 0.11916112 -0.138835389  0.3225370374
## 143 -0.3902964158 0.19273541 -0.760244820  0.0050620720
## 144 -0.0616769496 0.17640899 -0.399764444  0.3020565527
## 145  0.5779710850 0.13955182  0.302573718  0.8584200000
## 146  0.0280784726 0.17537969 -0.309256493  0.3858486609
## 147 -0.2371297993 0.18344687 -0.606346506  0.1244901972
## 148  0.0870711048 0.16800176 -0.238324337  0.4189615310
## 149 -0.0907731778 0.17126626 -0.431438071  0.2476504667
## 150  0.2106711034 0.16247436 -0.097624254  0.5286418315
## 151 -0.1004320154 0.18758487 -0.457926976  0.2777833242
## 152  0.2927863052 0.15888131 -0.018064912  0.6011312907
## 153 -0.1545042131 0.22650483 -0.606463637  0.2775133674
## 154 -0.0264660022 0.17300282 -0.356313864  0.3169644086
## 155  0.1301131398 0.22506491 -0.321798593  0.5721739949
## 156  0.1224206831 0.19227708 -0.249215358  0.5004368641
## 157 -0.0333716278 0.19171514 -0.412300803  0.3442530357
## 158 -0.4228169852 0.20043519 -0.800819845 -0.0191306201
## 159  0.2956023274 0.11152189  0.076255806  0.5079326310
## 160  0.1199851938 0.17083018 -0.209565566  0.4631624367
## 161  0.1904240641 0.19367429 -0.194696802  0.5750781849
## 162 -0.0206600602 0.23017397 -0.482701907  0.4277424230
## 163  0.1120462022 0.18061824 -0.240295068  0.4633206393
## 164 -0.2209541799 0.15978282 -0.528431876  0.1045518068
## 165 -0.2745546264 0.15657352 -0.574591088  0.0384209832
## 166 -0.1774740588 0.20539476 -0.573909889  0.2305886376
## 167 -0.3996230768 0.20802772 -0.815959167 -0.0019172073
## 168 -0.1195826315 0.19180422 -0.480267888  0.2561909117
## 169 -0.4101439050 0.20528598 -0.806985372  0.0005902036
## 170 -0.0620672062 0.21623006 -0.499300057  0.3623345796
## 171  0.1244439964 0.15516673 -0.177465135  0.4242820105
## 172 -0.2263308310 0.18624139 -0.582620760  0.1517117782
## 173 -0.1585357784 0.19511254 -0.534592772  0.2308130074
## 174  0.1968349933 0.13817580 -0.073441890  0.4650828446
## 175 -0.0482128702 0.17506467 -0.381965495  0.2991999676
## 176  0.0851508885 0.18338304 -0.286615235  0.4475140280
## 177 -0.0081811540 0.18765052 -0.370320553  0.3566915986
## 178 -0.0339306420 0.21926013 -0.443563950  0.3935362849
## 179  0.0642521768 0.17786358 -0.286624853  0.4189446593
## 180  0.1618294947 0.17748608 -0.180201132  0.5117088839
## 181  0.1625767129 0.18615122 -0.207551211  0.5315940904
## 182 -0.3479304861 0.19552153 -0.725241810  0.0401714921
## 183 -0.0736604057 0.17428417 -0.416273193  0.2772766827
## 184  0.2540727565 0.16680374 -0.058414707  0.5794173474
## 185 -0.4466473970 0.15431782 -0.743691984 -0.1447474564
## 186  0.1816633541 0.17358930 -0.161136683  0.5281742937
## 187 -0.0646187450 0.19276908 -0.448970544  0.3178217721
## 188  0.1181387258 0.21892963 -0.321171158  0.5407615704
## 189 -0.0922801212 0.18166460 -0.441773665  0.2658676452
## 190  0.0633306955 0.15387854 -0.233110516  0.3687586171
## 191 -0.1373840834 0.23723405 -0.609622218  0.3141880163
## 192 -0.1927435908 0.15317478 -0.491651075  0.1120770374
## 193  0.0811575039 0.11852846 -0.152096245  0.3158335355
## 194  0.0807796515 0.21216810 -0.334727711  0.4956162325
## 195  0.0509305092 0.17275981 -0.288786245  0.3921149623
## 196 -0.0491551314 0.20973221 -0.468257253  0.3682926729
## 197 -0.0623228342 0.17888812 -0.420455539  0.2940733538
## 198 -0.0492830855 0.17757984 -0.383707585  0.3059193449
## 199 -0.0206727593 0.19371728 -0.389861994  0.3623877728
## 200  0.0145170086 0.15508906 -0.293813276  0.3129378606
## 201  0.0796132012 0.18062972 -0.276553627  0.4362049287
## 202 -0.1045141305 0.19433842 -0.482472734  0.2829670497
## 203  0.0693181208 0.21960330 -0.371028090  0.4875462068
## 204 -0.0392898317 0.12420466 -0.274512196  0.2066093769
## 205  0.1441974819 0.16059943 -0.169800547  0.4596661226
## 206  0.1365505814 0.18531229 -0.230287008  0.4980182798
## 207 -0.1430139917 0.24211713 -0.636785110  0.3171161701
## 208  0.1239696191 0.14077466 -0.137095709  0.4148985120
## 209 -0.0404057662 0.15878268 -0.343553360  0.2716931833
## 210  0.0152965836 0.19974663 -0.378281070  0.4053682606
## 211  0.2370621149 0.19336321 -0.135915679  0.6293054614
## 212  0.1234657997 0.17688154 -0.237019425  0.4825573155
## 213  0.7993602339 0.13045864  0.545301149  1.0533841437
## 214 -0.2808050268 0.14216513 -0.551017524  0.0043331936
## 215  0.3970104453 0.16974052  0.076567809  0.7377096010
## 216  0.1596676758 0.20316334 -0.251198610  0.5500034244
## 217 -0.0385374160 0.17591268 -0.374317567  0.3152131888
## 218 -0.0648061543 0.12384106 -0.308116515  0.1814681150
## 219  0.3733816294 0.17052691  0.036981833  0.7115950373
## 220  0.1887070065 0.13072049 -0.068068576  0.4453845169
## 221  0.1080027206 0.16150839 -0.199360477  0.4228191472
## 222  0.3300621203 0.17189878 -0.003518501  0.6631077438
## 223 -0.0370251244 0.19280289 -0.417240784  0.3586792993
## 224  0.1215593725 0.21070079 -0.293535783  0.5298593989
## 225  0.0001376661 0.19092571 -0.374387998  0.3825161186
## 226  0.0114716994 0.17610308 -0.317790752  0.3604443533
## 227  0.0648224334 0.12748516 -0.171401152  0.3220237195
## 228  0.2530015766 0.16156383 -0.066889836  0.5698811854
## 229 -0.2245455627 0.18694014 -0.584784794  0.1463476499
## 230  0.0444833202 0.17569008 -0.296812387  0.4070484956
## 231  0.2075181109 0.17737098 -0.133724200  0.5694944911
## 232 -0.0319523623 0.18655161 -0.390764849  0.3366766343
## 233  0.0659900676 0.19887415 -0.325099851  0.4425550127
## 234  0.1604178197 0.17056050 -0.161509594  0.5033851350
## 235 -0.1051243524 0.23097465 -0.559794397  0.3535390386
## 236 -0.0795013269 0.20471109 -0.480381671  0.3249065780
## 237 -0.4245517079 0.20330685 -0.817885114 -0.0278367133
## 238  0.1304290727 0.17636568 -0.208951770  0.4758722799
## 239 -0.1490437518 0.17677423 -0.498109339  0.1942080804
## 240 -0.1110074823 0.18256591 -0.467384920  0.2613350103
## 241  0.2185181257 0.18810604 -0.145783529  0.5810614908
## 242  0.1195006813 0.19559490 -0.259318654  0.5139769684
## 243 -0.2371813839 0.15989859 -0.546214734  0.0714281398
## 244  0.1047311689 0.19977430 -0.298389424  0.4943182466
## 245 -0.0374733992 0.21854604 -0.481433792  0.3985435017
## 246  0.0771207920 0.19047265 -0.289275634  0.4484121918
## 247 -0.0122211566 0.18253186 -0.359257564  0.3565818365
## 248  0.0605685164 0.17802644 -0.278668625  0.4155522641
## 249 -0.0005168648 0.19368222 -0.372908868  0.3812901171
## 250  0.0286132601 0.21774012 -0.403866627  0.4454193123
## 251  0.0230901161 0.15075713 -0.264435627  0.3200914002
## 252 -0.1153306603 0.18760747 -0.493583158  0.2685406082
## 253  0.0157771712 0.17048066 -0.316343697  0.3475242233
## 254 -0.2127708258 0.21687482 -0.645321042  0.2036167407
## 255 -0.0012834033 0.17081699 -0.322090402  0.3478857662
## 256  0.1677722258 0.18346187 -0.191415533  0.5308891714
## 257 -0.0359045634 0.17491442 -0.366107100  0.3055155284
## 258 -0.0015081658 0.13240964 -0.245622627  0.2617007666
## 259  0.2801391732 0.18138632 -0.086876703  0.6477407355
## 260 -0.3774652251 0.22907431 -0.818310625  0.0690987116
## 261 -0.1755121121 0.18115389 -0.527552575  0.1927276291
## 262 -0.0514396781 0.14760907 -0.337317010  0.2398158762
## 263  0.0996025475 0.18813906 -0.277675417  0.4672962392
## 264  0.1077885206 0.16103479 -0.209602906  0.4210363940
## 265 -0.1097711740 0.16197516 -0.424064620  0.2058590228
## 266 -0.3021941785 0.18665708 -0.666166930  0.0645566649
## 267  0.2329867466 0.14873445 -0.050944018  0.5156686526
## 268 -0.0434172846 0.20765550 -0.451204177  0.3700951137
## 269  0.1512964253 0.14294859 -0.129580613  0.4382402589
## 270  0.0834044710 0.22472399 -0.362941898  0.5305248665
## 271  0.0685385192 0.19040056 -0.308649866  0.4397011746
## 272 -0.1888302347 0.20319486 -0.586665144  0.2070068888
## 273  0.3062175447 0.11344016  0.085578151  0.5300417333
## 274 -0.1271067828 0.19027078 -0.494201935  0.2407123937
## 275 -0.6011244025 0.21075555 -1.016501344 -0.1951673822
## 276 -0.0864816517 0.19266448 -0.471399050  0.2976468670
## 277  0.0327910046 0.18149021 -0.323917015  0.3872737822
## 278 -0.0262312216 0.22080181 -0.474313073  0.3984985193
## 279 -0.1076307054 0.18295255 -0.469864910  0.2460084865
## 280 -0.3229701888 0.19123177 -0.685278573  0.0492622132
## 281 -0.0780169427 0.19345025 -0.456631911  0.2962612719
## 282  0.2184604782 0.17754344 -0.130158453  0.5663374693
## 283 -0.0910956640 0.22303187 -0.529544020  0.3276333775
## 284 -0.0494987791 0.10231005 -0.252716733  0.1512267149
## 285 -0.1238376717 0.15649969 -0.423122365  0.1898349832
## 286  0.3217765779 0.18455851 -0.026096165  0.6873286691
## 287  0.0586445852 0.22380838 -0.382454505  0.5047923697
## 288  0.1073445077 0.19515052 -0.272727444  0.4887941101
## 289  0.0570183658 0.18738369 -0.306869169  0.4264430708
## 290  0.1124009366 0.21550298 -0.314284859  0.5330424696
## 291 -0.0637696180 0.20880190 -0.468289850  0.3481538023
## 292  0.1201383213 0.14773009 -0.164091435  0.4133232738
## 293  0.1149968423 0.19116115 -0.256679217  0.4871328245
## 294  0.0346421725 0.18024721 -0.302596268  0.4023576893
## 295  0.0220438483 0.17855853 -0.320213427  0.3751926057
## 296  0.2519958736 0.17609049 -0.082855223  0.5963938637
## 297 -0.1290022422 0.20287485 -0.521234800  0.2771617200
## 298 -0.0749869932 0.18638208 -0.436541164  0.3050707107
## 299  0.0423128790 0.12368617 -0.187854590  0.2917146830
## 300 -0.0779709669 0.20618174 -0.479516167  0.3173880510
## 301  0.4420967054 0.15463414  0.141901567  0.7373718469
## 302 -0.1578393090 0.22634280 -0.604383526  0.2951820208
## 303  0.0387421662 0.18952027 -0.313821169  0.4171134935
## 304 -0.4320638629 0.21491651 -0.861925230 -0.0072848971
## 305 -0.1889355704 0.14033079 -0.457763724  0.0947893996
## 306 -0.3963274913 0.22844804 -0.845554936  0.0433969668
## 307 -0.1054267715 0.17838984 -0.450865812  0.2513355034
## 308 -0.0001082094 0.20524215 -0.406800966  0.4040224741
## 309 -0.0595438907 0.21253913 -0.479072485  0.3573142088
## 310 -0.1573862324 0.18486885 -0.508342754  0.2270454079
## 311 -0.0414743001 0.18132168 -0.388951920  0.3241156581
## 312  0.0014058383 0.14259846 -0.264617667  0.2887558047
## 313  0.1052024006 0.20190611 -0.285545619  0.5081615651
## 314  0.0903725761 0.18528892 -0.276806580  0.4619000346
## 315  0.2776046419 0.14733363 -0.011409904  0.5697374765
## 316 -0.0866789986 0.23283294 -0.548629806  0.3587074705
## 317  0.0588554730 0.21670742 -0.384247636  0.4899515600
## 318  0.1164732474 0.14915296 -0.166325683  0.4044292669
## 319 -0.1743558445 0.14327513 -0.451192669  0.1187274472
## 320 -0.1344955350 0.19044118 -0.505882903  0.2416399085
## 321  0.2389947492 0.16443739 -0.072962985  0.5602739157
## 322 -0.0799039107 0.23459904 -0.534676073  0.3716840787
## 323  0.1476526691 0.18785588 -0.223723136  0.5137101617
## 324  0.4484865212 0.13135676  0.188776629  0.7056161654
## 325 -0.1566584992 0.18468461 -0.504839240  0.2117585145
## 326  0.2687602903 0.20725491 -0.125698502  0.6687939373
## 327 -0.0814619674 0.17768514 -0.420528145  0.2739747476
## 328  0.0760654049 0.10747632 -0.138463549  0.2930488374
## 329 -0.1342247156 0.15130924 -0.435793088  0.1587705220
## 330 -0.0728407358 0.18653895 -0.430381670  0.2911833220
## 331 -0.1161370906 0.18091461 -0.462795006  0.2410515166
## 332  0.0583698270 0.17357124 -0.281925437  0.4016604086
## 333 -0.0517663538 0.19672595 -0.435341104  0.3305421136
## 334  0.1136251007 0.18908401 -0.258233664  0.4903655790
## 335  0.1510916462 0.14360619 -0.131584436  0.4302893534
## 336  0.1717248013 0.18216338 -0.176665683  0.5405938881
## 337 -0.4628162536 0.19574086 -0.846184966 -0.0783695471
## 338 -0.0748698419 0.18195983 -0.436464614  0.2881943324
## 339 -0.1876709301 0.17794830 -0.539320493  0.1743914572
## 340 -0.2579774723 0.20432331 -0.649025502  0.1418369443
## 341  0.2104145989 0.18457034 -0.143321224  0.5739650514
## 342  0.1776001363 0.16989483 -0.152753060  0.5189153337
## 343  0.0654974122 0.17925484 -0.274248180  0.4196219570
## 344 -0.1017169582 0.18083659 -0.441694898  0.2523512349
## 345  0.1048394955 0.17014964 -0.218809334  0.4509037323
## 346 -0.0325456923 0.18188673 -0.378803059  0.3286011409
## 347 -0.0121342218 0.19362359 -0.381628816  0.3809333574
## 348  0.0167630209 0.18978636 -0.349020391  0.4135472431
## 349  0.1483694998 0.12176403 -0.087207990  0.3937840095
## 350  0.0205807532 0.17069578 -0.307100106  0.3585503167
## 351 -0.1338867245 0.15217490 -0.423902751  0.1666650915
## 352 -0.2538556965 0.16183031 -0.567380249  0.0761528008
## 353  0.2143970682 0.17663972 -0.126524536  0.5738078869
## 354  0.3272086843 0.15846169  0.019119693  0.6371132279
## 355 -0.3535274972 0.13630621 -0.615008667 -0.0854695332
## 356  0.0512159872 0.18191359 -0.294236136  0.4021769597
## 357  0.3225783790 0.13153599  0.069037761  0.5872702868
## 358 -0.0038721777 0.15810662 -0.313482864  0.3168119280
## 359  0.2591472638 0.17308388 -0.079289365  0.6071790072
## 360 -0.2194700900 0.19032658 -0.581956486  0.1560298392
## 361  0.3578741886 0.16608660  0.047791055  0.6883511908
## 362  0.0269754081 0.18982652 -0.347389204  0.3989572806
## 363 -0.0159705501 0.17785608 -0.361048767  0.3321765315
## 364 -0.3000818654 0.14544546 -0.585756329 -0.0111342600
## 365  0.1259849841 0.17064291 -0.204868965  0.4649082593
## 366  0.0421044614 0.17366272 -0.290402165  0.3876961280
## 367  0.0414415434 0.21686024 -0.384147520  0.4522906002
## 368  0.0271058727 0.14525158 -0.255579421  0.3109970977
## 369 -0.1341103755 0.19491990 -0.505835463  0.2526299466
## 370 -0.2366723519 0.20850434 -0.648665124  0.1629715372
## 371  0.1352929652 0.20362984 -0.267942505  0.5328948935
## 372 -0.0854516488 0.20027541 -0.482247950  0.3058128919
## 373 -0.0169358445 0.18909692 -0.389924283  0.3456643458
## 374 -0.0273530894 0.15874714 -0.336027325  0.2974528051
## 375  0.0686170018 0.17430081 -0.263450839  0.4107595343
## 376 -0.0317809377 0.15936954 -0.352255680  0.2801606310
## 377  0.2497669846 0.15685899 -0.057882392  0.5607368088
## 378 -0.1764011226 0.17577661 -0.510374198  0.1818608100
## 379 -0.2743464949 0.15582482 -0.574639474  0.0286385727
## 380 -0.0706868400 0.21688054 -0.496066631  0.3505679590
## 381  0.0564726902 0.14464252 -0.221718662  0.3386276825
## 382  0.0505958521 0.16553564 -0.268566754  0.3753622454
## 383 -0.1964763785 0.17771862 -0.541909389  0.1584338136
## 384 -0.0554189466 0.12590426 -0.298214656  0.1908490038
## 385 -0.0805826415 0.11664032 -0.311130649  0.1559270953
## 386  0.2196120483 0.17650616 -0.118111886  0.5626378851
## 387 -0.1671511773 0.19262219 -0.542201196  0.2271401032
## 388 -0.2297134623 0.19999921 -0.628470694  0.1674988838
## 389  0.1045250501 0.17937709 -0.246463588  0.4718738787
## 390 -0.4299439827 0.23410162 -0.888677468  0.0236958680
## 391 -0.0116400560 0.19750100 -0.388330711  0.3674632036
## 392 -0.2170591033 0.18085100 -0.572080153  0.1517133142
## 393 -0.3271456595 0.12319084 -0.561048762 -0.0774856637
## 394 -0.0247353057 0.17811061 -0.372257269  0.3381001393
## 395  0.0648330886 0.10394298 -0.139959714  0.2696394442
## 396  0.5361750830 0.14675137  0.248309648  0.8228407472
## 397  0.0482278849 0.17652474 -0.297650554  0.3961897740
## 398  0.0554019177 0.17589469 -0.283341103  0.3954507710
## 399 -0.0578397523 0.15700691 -0.356509941  0.2493616008
## 400 -0.2179237312 0.19688216 -0.600175893  0.1667226296
## 401 -0.0629972577 0.18558145 -0.425787911  0.3032173399
## 402 -0.1655554910 0.22594962 -0.606377310  0.2687672232
## 403 -0.1802020652 0.14331094 -0.459153897  0.1098065334
## 404  0.1798866284 0.21734974 -0.237715582  0.5991797666
## 405 -0.2655428258 0.11025845 -0.479713665 -0.0463348139
## 406  0.2902957935 0.16965179 -0.036367868  0.6294289382
## 407  0.0360143828 0.18122800 -0.318537712  0.3982451013
## 408  0.5395438187 0.11589425  0.314647338  0.7636119101
## 409  0.3352232338 0.20306217 -0.077855644  0.7222663614
## 410  0.2383670912 0.09431584  0.048161955  0.4205543965
## 411 -0.2302644241 0.17551719 -0.571028453  0.1220599080
## 412 -0.0403090045 0.17467285 -0.374978047  0.3110468160
## 413 -0.0290890933 0.18003185 -0.380252188  0.3354353854
## 414 -0.0328681143 0.17229882 -0.362120782  0.3070958398
## 415 -0.1721931508 0.15157637 -0.467951265  0.1263732964
## 416 -0.0825107229 0.18039063 -0.437768923  0.2788616467
## 417  0.1039749180 0.17490457 -0.235955676  0.4531391858
## 418  0.4375115657 0.16903151  0.104517020  0.7780547614
## 419 -0.0533376867 0.22886737 -0.503868665  0.3831130621
## 420  0.2247468219 0.18800099 -0.143058882  0.6024775413
## 421  0.0145643218 0.17566428 -0.321330954  0.3735308791
## 422 -0.3653218876 0.19607971 -0.745389321  0.0184948522
## 423 -0.0270593400 0.20990606 -0.433456788  0.3859864586
## 424  0.0936675844 0.18171667 -0.254571402  0.4538820168
## 425  0.0249990179 0.16334059 -0.288775998  0.3537644803
## 426 -0.1217078061 0.12389557 -0.363950769  0.1222006530
## 427 -0.3926852785 0.19493303 -0.769500425 -0.0107900630
## 428  0.3727253522 0.20567819 -0.034291323  0.7726814376
## 429  0.0834539388 0.22461595 -0.366238941  0.5250481617
## 430  0.1394749863 0.16171876 -0.175592701  0.4727368900
## 431 -0.1425180894 0.17006074 -0.480760139  0.1860089843
## 432  0.0946487485 0.21324505 -0.329111167  0.5135655239
## 433  0.3316563519 0.18770301 -0.041190832  0.6999919363
## 434  0.0757694428 0.21932699 -0.357096850  0.4931382957
## 435 -0.0182238461 0.17195500 -0.354607496  0.3216457000
## 436 -0.0769020685 0.18513371 -0.426108079  0.2921889376
## 437 -0.1755314355 0.20289185 -0.571985725  0.2161480979
## 438  0.1324274456 0.13314069 -0.124812961  0.4006615330
## 439  0.3657999913 0.15130228  0.075084200  0.6649742539
## 440  0.0889671475 0.15369904 -0.206358347  0.3843263980
## 441  0.2428054065 0.17502993 -0.097049219  0.5903521613
## 442  0.2822716029 0.18614535 -0.076206774  0.6550703281
## 443  0.0346905669 0.15230566 -0.263149528  0.3285542588
## 444 -0.0609217973 0.20206736 -0.449761484  0.3380484813
## 445 -0.1584642048 0.18310847 -0.506322088  0.2012590952
## 446  0.0685544458 0.21285760 -0.344788902  0.4848387066
## 447 -0.0132767485 0.10459904 -0.218023247  0.1934261038
## 448 -0.1349009158 0.19343644 -0.519980522  0.2540439468
## 449 -0.2598355967 0.14215878 -0.540434849  0.0174749179
## 450  0.0299300313 0.13100582 -0.220177727  0.2888079576
## 451 -0.0086154365 0.21139196 -0.422248537  0.3992881881
## 452  0.3369762146 0.16995936  0.009841266  0.6745081965
## 453 -0.3479477786 0.18088866 -0.699942078  0.0019600296
## 454  0.1569696898 0.11959115 -0.082430083  0.3888728923
## 455  0.1509617219 0.16369325 -0.164342496  0.4775480561
## 456 -0.0148926129 0.18163046 -0.379105037  0.3416312758
## 457  0.1662208743 0.20755237 -0.237379802  0.5736405422
## 458  0.0906669877 0.17358237 -0.243985207  0.4241801630
## 459 -0.2561742376 0.19115690 -0.629036590  0.1270673737
## 460  0.2467217420 0.16940140 -0.084396228  0.5915892858
## 461  0.2035625315 0.13930376 -0.062811518  0.4754905728
## 462  0.0581428707 0.15786305 -0.256202791  0.3709647704
## 463  0.0092967320 0.17696913 -0.326062253  0.3695007524
## 464  0.1173167564 0.20367478 -0.266978371  0.5139405320
## 465 -0.0551547682 0.23458727 -0.530645534  0.4103237713
## 466  0.0335384007 0.15705180 -0.263548580  0.3459392062
## 467  0.1293466194 0.14171244 -0.147084274  0.4050332867
## 468 -0.0281452702 0.15214174 -0.325285289  0.2719680251
## 469  0.0233193958 0.16771330 -0.300066166  0.3442837980
## 470 -0.0625439023 0.19150655 -0.429549142  0.3273313173
## 471  0.1518110106 0.13646647 -0.107124036  0.4259576907
## 472  0.3744774571 0.11108608  0.157263273  0.5906163539
## 473  0.2666050448 0.12965173  0.023440781  0.5308201966
## 474 -0.4471897560 0.13402481 -0.709859884 -0.1854968400
## 475  0.1372635270 0.18828645 -0.233775147  0.5059411152
## 476  0.4886379894 0.17210224  0.154647544  0.8350165264
## 477 -0.0235101736 0.10899707 -0.238020270  0.1875242045
## 478  0.0828956303 0.17354270 -0.260957652  0.4331562318
## 479  0.0927451599 0.14031038 -0.174629048  0.3733044694
## 480 -0.1930354085 0.17594504 -0.531180285  0.1528401995
## 481 -0.0671281851 0.18984763 -0.443736535  0.3062769712
## 482  0.2282475535 0.10572059  0.027157252  0.4348654229
## 483  0.1713792088 0.09701267 -0.016523798  0.3613811600
## 484 -0.0060980989 0.13232052 -0.262396659  0.2585625272
## 485 -0.0058269305 0.13419179 -0.269737863  0.2640285677
## 486  0.0895035051 0.20178748 -0.301519472  0.4879072437
## 487  0.1735001086 0.14354851 -0.102052851  0.4679748618
## 488 -0.1349133984 0.12703433 -0.380045361  0.1149490633
## 489  0.3518564073 0.13511313  0.091571496  0.6240924458
## 490 -0.0180988403 0.19325178 -0.392937297  0.3628717343
## 491  0.0415636966 0.20704870 -0.356678935  0.4512493489
## 492 -0.3260655415 0.16287545 -0.642567127 -0.0015496434
## 493  0.2124176451 0.17265606 -0.118082918  0.5491537104
## 494  0.1208766157 0.16014767 -0.199698234  0.4415365841
## 495 -0.1458257341 0.16495877 -0.460947509  0.1737555893
## 496 -0.2181670144 0.18218566 -0.565551808  0.1313147040
## 497 -0.0661379536 0.16076296 -0.381791165  0.2425040342
## 498 -0.2043026854 0.23122290 -0.677616293  0.2489682763
## 499 -0.0090876737 0.17409915 -0.349676117  0.3437751411
## 500 -0.2191931244 0.23112077 -0.672030561  0.2290996598
## 501  0.1576671369 0.16996717 -0.174770266  0.4955281922
## 502  0.0563764186 0.18254043 -0.303073284  0.4315640237
## 503 -0.0290334036 0.16110512 -0.338796408  0.2904391254
## 504 -0.1242576913 0.13625984 -0.389834675  0.1411277994
## 505  0.2495602021 0.18731101 -0.111086428  0.6238145388
## 506  0.2682722053 0.17810214 -0.080178006  0.6139031843
## 507  0.0986297415 0.21759530 -0.325697119  0.5279442275
## 508  0.1548417351 0.16385439 -0.160026730  0.4747756039
## 509  0.1430618829 0.18765872 -0.223395603  0.5060628694
## 510  0.0228460993 0.20516296 -0.370096314  0.4230841325
## 511  0.1941077438 0.11514949 -0.027088478  0.4242489119
## 512  0.0483233406 0.10647226 -0.164529008  0.2533039768
## 513 -0.0690843631 0.17994580 -0.412440361  0.2861175167
## 514  0.0528468105 0.10640239 -0.152209566  0.2625942793
## 515  0.0945599473 0.17957950 -0.254885802  0.4525856986
## 516  0.1559683280 0.18113943 -0.190556013  0.5123562927
## 517 -0.0600384892 0.19552699 -0.442905170  0.3242395593
## 518  0.1064225425 0.22010752 -0.319475285  0.5326075437
## 519 -0.2064662512 0.14849359 -0.493295216  0.0860526166
## 520  0.1323194067 0.12760958 -0.118045989  0.3856363289
## 521 -0.3025931902 0.17483373 -0.643289927  0.0419292555
## 522 -0.1391959722 0.15099162 -0.428127073  0.1660905339
## 523 -0.0028076710 0.13869530 -0.268329527  0.2684357993
## 524  0.1558579733 0.19491292 -0.218140734  0.5280538759
## 525  0.0313790270 0.17620526 -0.307453080  0.3760927007
## 526 -0.0057257278 0.18535780 -0.364131318  0.3582760853
## 527  0.0256672619 0.14307689 -0.250370815  0.3117058667
## 528 -0.0658071034 0.15273427 -0.368008830  0.2360065600
## 529 -0.4592205310 0.19523176 -0.831228715 -0.0756716857
## 530  0.0672525604 0.18133525 -0.296552228  0.4307296254
## 531  0.0781431860 0.17471671 -0.266217181  0.4154327894
## 532  0.1384962536 0.17586703 -0.199327139  0.4883255912
## 533  0.1994880591 0.21191508 -0.227092708  0.6284233936
## 534  0.0710924943 0.15697587 -0.238652310  0.3796024087
## 535 -0.3196858360 0.14854641 -0.609321972 -0.0278931749
## 536 -0.0632134438 0.14954356 -0.349239147  0.2353549041
## 537 -0.0609809940 0.18399428 -0.417700774  0.3037997360
## 538  0.2128043150 0.19333172 -0.165443242  0.5933424430
## 539 -0.1660645524 0.15485165 -0.464312551  0.1428510884
## 540 -0.1528013954 0.18325274 -0.508528127  0.2036169412
## 541  0.1265072141 0.16998650 -0.200192052  0.4677749444
## 542  0.0699246200 0.15173658 -0.228416207  0.3648176599
## 543  0.0732427063 0.17461565 -0.257298726  0.4094764501
## 544 -0.1666229584 0.21484954 -0.590487779  0.2546309873
## 545  0.1269309711 0.10740425 -0.078572075  0.3421223976
## 546 -0.0502690186 0.23278057 -0.515733273  0.4098306680
## 547 -0.2662026168 0.18501936 -0.626154520  0.1093763247
## 548  0.1540156944 0.17666703 -0.187293248  0.5065229694
## 549  0.0727423344 0.22805181 -0.377910832  0.5069716094
## 550 -0.4033143941 0.17923749 -0.747834713 -0.0606921637
## 551 -0.0614309343 0.18246818 -0.413696162  0.2932290588
## 552 -0.1197239946 0.18847480 -0.488943688  0.2460339636
## 553 -0.0525441101 0.17758162 -0.391252831  0.2983895506
## 554 -0.3869048975 0.19961297 -0.782523556 -0.0011821440
## 555 -0.1113050354 0.15248328 -0.402542509  0.1918725671
## 556 -0.1704309931 0.15727920 -0.485405233  0.1393523457
## 557  0.2566014706 0.13821437 -0.011923047  0.5299353928
## 558 -0.0115342421 0.18195575 -0.362718358  0.3557450899
## 559 -0.0142792914 0.14626634 -0.292924769  0.2829644569
## 560  0.0165770701 0.19737793 -0.371195690  0.4011741456
## 561 -0.1811082046 0.14864515 -0.458357120  0.1096142783
## 562 -0.0584432417 0.17974431 -0.399524600  0.3075509859
## 563  0.0527858782 0.17357353 -0.284756216  0.3971592981
## 564  0.2108832029 0.13367943 -0.054420559  0.4713030093
## 565  0.1881459809 0.18764347 -0.166332675  0.5601477043
## 566 -0.0923294029 0.17696004 -0.430528677  0.2703231058
## 567 -0.4811990761 0.17762166 -0.824925343 -0.1359414960
## 568 -0.1673693638 0.16879849 -0.493472747  0.1660415041
## 569 -0.2262493966 0.20338587 -0.618881362  0.1662996032
## 570 -0.1727459677 0.13305745 -0.431557126  0.0908122651
## 571  0.2175315896 0.17069036 -0.102123005  0.5560164725
## 572 -0.0569763032 0.20996496 -0.468227300  0.3464440695
## 573 -0.1255988914 0.17870701 -0.475025152  0.2356081917
## 574  0.2267787925 0.19345407 -0.161377650  0.6014291423
## 575 -0.2646790021 0.20037662 -0.655793956  0.1278281694
## 576  0.0751011506 0.17160613 -0.260784658  0.4170682763
## 577  0.2027213651 0.17787780 -0.142192474  0.5514624107
## 578  0.2786002602 0.17123350 -0.061436916  0.6171669659
## 579 -0.0363448273 0.20825149 -0.440901532  0.3756865337
## 580 -0.0732286458 0.18185174 -0.421318949  0.2882629890
## 581  0.3189821084 0.16840063 -0.001615695  0.6601210763
## 582 -0.2306614991 0.19681921 -0.610484352  0.1585452332
## 583 -0.2769049443 0.15691900 -0.582731419  0.0346119327
## 584 -0.5306507320 0.21581152 -0.951415010 -0.0940125180
## 585  0.0265795689 0.18452320 -0.325954685  0.3930073984
## 586 -0.1988428157 0.16549785 -0.524172209  0.1343016793
## 587  0.1936056628 0.19129940 -0.169984994  0.5645758228
## 588  0.0789512968 0.18403852 -0.274657282  0.4389615787
## 589 -0.1103735074 0.19070159 -0.479976884  0.2691591917
## 590  0.0390055045 0.16058800 -0.266409857  0.3559416349
## 591  0.0457944987 0.20668418 -0.350428682  0.4458098244
## 592  0.1497100344 0.13793624 -0.106539902  0.4224347166
## 593 -0.1185879461 0.16004502 -0.424038832  0.2055308580
## 594  0.0576597293 0.13725973 -0.205838226  0.3319339839
## 595 -0.0328119643 0.18182601 -0.384495129  0.3327962110
## 596 -0.1127241689 0.18237265 -0.466572617  0.2501785452
## 597  0.0396077852 0.17265921 -0.286481918  0.3860878442
## 598  0.0621885964 0.16914812 -0.273397393  0.3940801773
## 599 -0.3009054095 0.21340176 -0.711882622  0.1176436023
## 600  0.5041118040 0.16120959  0.191316503  0.8155616797
## 601  0.2508658155 0.17631567 -0.088584605  0.6045509120
## 602  0.5585660345 0.15283528  0.262077639  0.8521446285
## 603 -0.3398272822 0.16714605 -0.658841172 -0.0186131220
## 604 -0.0035377483 0.18134648 -0.350900205  0.3684636723
## 605  0.1164782569 0.18191233 -0.232176815  0.4809791060
## 606 -0.1188092620 0.18387549 -0.472476341  0.2555177623
## 607 -0.0121520517 0.17217711 -0.346374572  0.3348067832
## 608  0.0274542034 0.18477250 -0.330523116  0.4073352238
## 609 -0.0717046632 0.20199743 -0.465443378  0.3342715610
## 610  0.0686527255 0.16667979 -0.245812062  0.4000065034
## 611  0.0146392916 0.20610040 -0.392514378  0.4184032895
## 612 -0.0974465760 0.15013175 -0.383651498  0.2066511780
## 613 -0.2414034439 0.22487993 -0.681250558  0.1953352885
## 614  0.0422872473 0.12713128 -0.202725045  0.2923445656
## 615  0.0366311532 0.16984690 -0.292388536  0.3802022151
## 616  0.0194550849 0.13419867 -0.240813614  0.2889459045
## 617  0.0974529330 0.17900761 -0.251888257  0.4520738054
## 618  0.0095281607 0.17620896 -0.335874540  0.3560806626
## 619  0.3318479154 0.12829296  0.086772739  0.5849586846
## 620  0.1157689807 0.17022212 -0.223334684  0.4533637024
## 621 -0.2495778355 0.18669957 -0.618049347  0.1252135361
## 622 -0.0891369279 0.23443644 -0.564212902  0.3685648422
## 623  0.1426922154 0.18421972 -0.208736856  0.5040072462
## 624  0.2092297278 0.17517355 -0.123712442  0.5521495064
## 625  0.1007421558 0.18433925 -0.246777808  0.4629565796
## 626 -0.0190000463 0.22795112 -0.459524105  0.4332747884
## 627  0.0218129187 0.16803095 -0.306854208  0.3519023023
## 628 -0.2800901424 0.19748351 -0.671095914  0.0988902558
## 629 -0.0962001568 0.13278165 -0.356418159  0.1797605640
## 630  0.0568385129 0.19131325 -0.324907392  0.4489713924
## 631  0.2422308281 0.22119809 -0.198802252  0.6809657642
## 632  0.1434334275 0.13166424 -0.110173967  0.4050487559
## 633  0.0341462350 0.17971270 -0.316344267  0.3863658761
## 634  0.0383724114 0.18929867 -0.321389712  0.4081233509
## 635 -0.1679583713 0.15901920 -0.477457651  0.1495383879
## 636  0.0800080332 0.17710478 -0.254625291  0.4298422985
## 637  0.1125379424 0.20460572 -0.293970471  0.5097337642
## 638 -0.3668789717 0.19394766 -0.749159792  0.0187820158
## 639 -0.1335647874 0.20000020 -0.524834738  0.2462672988
## 640 -0.0133483138 0.18840948 -0.374869371  0.3613036004
## 641  0.1485627022 0.17415320 -0.183534883  0.4863550470
## 642 -0.1326152589 0.22646423 -0.562183781  0.3089813780
## 643 -0.0454823183 0.16157993 -0.356630398  0.2760858587
## 644 -0.1913853791 0.18534575 -0.550573088  0.1791233588
## 645  0.0348765877 0.19242796 -0.338961940  0.4255241029
## 646 -0.2775772675 0.15082268 -0.563485173  0.0178555194
## 647 -0.1791153494 0.17285562 -0.517479540  0.1673827616
## 648  0.0936837970 0.11229035 -0.129158227  0.3097761145
## 649 -0.1187097855 0.19769469 -0.501385094  0.2814722997
## 650 -0.0380658495 0.17858328 -0.383811545  0.3223290997
## 651 -0.1020070166 0.13365889 -0.357676507  0.1536277754
## 652  0.2364805285 0.19452494 -0.128024940  0.6141933096
## 653  0.2028795589 0.17139172 -0.131494940  0.5441813231
## 654 -0.3825264970 0.14361331 -0.666116254 -0.0963726986
## 655  0.1543662324 0.18585214 -0.203520133  0.5238198805
## 656 -0.2414987622 0.20154199 -0.639144419  0.1584920925
## 657 -0.0664547503 0.17344532 -0.399360861  0.2701181021
## 658  0.3219509003 0.15993250  0.013230391  0.6330765957
## 659 -0.3265902436 0.20242950 -0.725570200  0.0649020414
## 660 -0.6241963220 0.15816923 -0.920752220 -0.3012322285
## 661  0.1035755188 0.15363172 -0.200649808  0.4060093105
## 662  0.5401033030 0.14469926  0.261308000  0.8301252004
## 663  0.0989854031 0.16626630 -0.232454673  0.4285512995
## 664  0.0551983944 0.17316877 -0.269787114  0.4084912775
## 665 -0.0904449147 0.16468037 -0.408373853  0.2393601013
## 666 -0.3692944505 0.20453168 -0.783896751  0.0218171516
## 667 -0.0531822643 0.19667934 -0.439597268  0.3379446672
## 668  0.3062003633 0.14515834  0.017209151  0.5841723620
## 669  0.0146557443 0.22364506 -0.427793345  0.4498959091
## 670 -0.1686634647 0.20354168 -0.561731793  0.2316270373
## 671  0.0628846740 0.14655332 -0.213538865  0.3535806137
## 672 -0.0575052108 0.17840899 -0.407250428  0.2892070313
## 673 -0.2654374387 0.19700278 -0.647304025  0.1201401292
## 674  0.1299376327 0.19695018 -0.259240057  0.5257746373
## 675 -0.1847163721 0.19612263 -0.569052510  0.1995462123
## 676  0.2797740645 0.19398128 -0.095857261  0.6540761527
## 677 -0.0051186582 0.20669980 -0.412135056  0.4063285922
## 678 -0.1200784649 0.18954356 -0.480900963  0.2559597082
## 679  0.0501430461 0.18993015 -0.314596184  0.4346839874
## 680 -0.0779294320 0.17587742 -0.424907673  0.2802776735
## 681  0.2375359920 0.14534529 -0.045013754  0.5321092279
## 682 -0.1665470774 0.18585546 -0.537038856  0.2002514812
## 683  0.4125253559 0.13519675  0.154322546  0.6775874797
## 684 -0.1619417200 0.18026365 -0.508374236  0.2075669470
## 685 -0.1693819711 0.18771959 -0.537170475  0.2090437408
## 686 -0.0139283164 0.21287288 -0.426252528  0.3913983814
## 687 -0.0084184624 0.23269488 -0.465517761  0.4469265143
## 688  0.1191901158 0.17853257 -0.218713339  0.4796850794
## 689  0.2084526313 0.20113579 -0.181661613  0.6093063277
## 690 -0.2356191748 0.19811187 -0.623851288  0.1417272346
## 691  0.0272891665 0.19480450 -0.350155851  0.4078414751
## 692 -0.0838773583 0.14874369 -0.377078891  0.2124520953
## 693  0.0539351533 0.20872996 -0.350846006  0.4701650652
## 694  0.0962407036 0.19294716 -0.289392922  0.4616081015
## 695  0.1983792169 0.17017430 -0.128124533  0.5413238555
## 696 -0.2111931062 0.18598920 -0.561891105  0.1634088511
## 697  0.0237657661 0.19122123 -0.352302967  0.3944676782
## 698 -0.1819055952 0.19664344 -0.573170633  0.2033567951
## 699 -0.0052401170 0.18619301 -0.360136728  0.3695409981
## 700  0.1118987388 0.17029356 -0.221662881  0.4548703377
## 701  0.0123396555 0.21738355 -0.411909300  0.4330043520
## 702  0.0391676777 0.17978220 -0.312409463  0.3985103351
## 703 -0.1217359693 0.20274636 -0.519597879  0.2865060757
## 704  0.0105157624 0.18474716 -0.348140450  0.3802750164
## 705  0.3434247412 0.13680712  0.079256831  0.6264496239
## 706 -0.2198581904 0.23337020 -0.682391348  0.2345889803
## 707  0.4410793144 0.16737753  0.114249171  0.7652874854
## 708 -0.3053350625 0.21464366 -0.737180921  0.1276122222
## 709 -0.1770030110 0.23475317 -0.649072006  0.2812344791
## 710 -0.2423635931 0.20965216 -0.669925535  0.1666215829
## 711  0.3048248240 0.18563136 -0.048186671  0.6720541707
## 712  0.1805467420 0.11080933 -0.029981576  0.3999608315
## 713  0.1343806783 0.17604926 -0.204706593  0.4929618936
## 714  0.2266350778 0.14704319 -0.055136811  0.5303081041
## 715  0.2790128111 0.16846154 -0.043562529  0.6136865456
## 716  0.0473530236 0.15379595 -0.249689996  0.3511284412
## 717  0.0652903784 0.16068525 -0.248353408  0.3788352574
## 718 -0.1320200536 0.18054626 -0.484495441  0.2317348994
## 719 -0.1734654838 0.17078191 -0.498227143  0.1663123221
## 720  0.2890695116 0.13802853  0.024823923  0.5606640402
## 721  0.0810943203 0.14542566 -0.207494695  0.3656125250
## 722 -0.0690594402 0.17815570 -0.406727703  0.2866168726
## 723  0.1375596003 0.19010665 -0.227173389  0.4976921464
## 724 -0.0511941292 0.15885418 -0.357653472  0.2551134396
## 725  0.1638724021 0.15508584 -0.135592730  0.4731671847
## 726 -0.0463028384 0.18171111 -0.404932931  0.3121906502
## 727  1.0142589951 0.13270125  0.761311029  1.2794111823
## 728 -0.1890575857 0.15713081 -0.497767204  0.1230104006
## 729 -0.3022445414 0.20076295 -0.689567118  0.0950660967

References

Dvorzak, M., and H. Wagner. 2016. “Sparse Bayesian Modelling of Underreported Count Data.” Statistical Modelling 16 (1): 24–46. https://doi.org/10.1177/1471082x15588398.

Paananen, Topi, Juho Piironen, Paul-Christian Bürkner, and Aki Vehtari. 2019. “Implicitly Adaptive Importance Sampling.” arXiv e-Prints, June, arXiv:1906.08850. http://arxiv.org/abs/1906.08850.

Ray, B., D. Posnett, V. Filkov, and P. T. Devanbu. 2014. “A Large Scale Study of Programming Languages and Code Quality in GitHub.” In SIGSOFT Symposium on the Foundations of Software Engineering, 155–65. ACM.

Vehtari, A., A. Gelman, and J. Gabry. 2017. “Practical Bayesian Model Evaluation Using Leave-One-Out Cross-Validation and WAIC.” Statistics and Computing 27: 1413–32. https://doi.org/10.1007/s11222-016-9696-4.

---

1. Wikipedia entry for the negative binomial distribution↩︎

2. https://mc-stan.org/docs/2_23/reference-manual/divergent-transitions.html↩︎

3. In this replication package we follow recommendations on priors. For an introduction please see (here)[https://github.com/stan-dev/stan/wiki/Prior-Choice-Recommendations].↩︎

4. We used Hamiltonian Monte Carlo with the No U-Turn Sampler (NUTS), which brms uses thanks to the functionality from the Stan package.↩︎

5. In the remaining executions we’ll refrain from showing the output, but we have made sure that no warnings exist when sampling (simply running this Rmd file yourself will allow you to check this).↩︎

6. In the remainder of the replication package we will not mention this anymore, but this should always be done before we put any trust in a model.↩︎

7. A very good introduction to LKJ can be read here: http://srmart.in/is-the-lkj1-prior-uniform-yes/↩︎

8. For an excellent introduction to these terms please see this post.↩︎

9. Credible intervals↩︎

10. Projection predictive variable selection↩︎